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Use of Continued Fractions in High 
Speed Computing 


1. Introduction. In the course of carrying out computations required 
for numerical solution of problems it is frequently necessary to have avail- 
able the value of one or more functions for various values of their arguments. 
If a high speed machine is being used, it is usually not efficient to look up 
tables of the functions outside the machine or store them in the internal 
memory and it is therefore necessary to calculate the values of the function 
whenever they are required. This is generally done either by using a rational 
approximation to the function or a finite number of terms of an infinite 
process. Usually the infinite process used is that of the power series. The 
purpose of this paper is to show that another infinite process, a continued 
fraction expansion, may, in some cases, be more efficient. 

The choice of which infinite process to use may depend on: 


1. Properties of the method of computing, such as 
(a) the number of orders required to program the calculation 
(b) rounding error involved 
(c) the magnitude of the numbers entering into the calculations 
2. Properties of the series, such as 
(d) speed of convergence 
(e) region of convergence 
3. Properties of the machine being used, such as 


(f) amount of high speed storage available 
(g) what factors determine the total time spent on a problem. 


Since methods of computing continued fractions are not as well known 
as they deserve to be, the following section will be devoted to examples of 
three methods of computation. In section 3 the speed of convergence of 
some particular continued fractions will be compared with that of the corre- 
sponding power series and in section 4 we shall discuss some situations in 
which continued fractions might be used. 

2. Computing Continued Fractions. A continued fraction expansion of 
the function f(x) takes the form 


(1) fix) = bo +a 

bb + a2 

be + <>, 
where the a’s and b’s may be functions of x. For convenience the expansion 
is usually written in the form 
(2) f(x) =bo +a a as 
by + be + bs + +>. 
127 
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The following are examples of such expansions 
(3) e = I zs £ £ & & 
i—-1+2-3+2-5+--- 
(4) Inx=x-1 P(x — 1) 12(x — 1) 2(x — 1) 
.. + 2 + 3 + 4 
Ae— 1) He —-1) 
+ 5 + 6 + .* 
) arctanx =x (x)? (2x) (3x) 
i+ 2+ & + F +:-- 
In practice a finite number of terms are used to approximate the function; 
the so-called “nth approximant” is given by 
(6) fr(x) =bo +a Gn 
bh +be+--- +d, 
Three different ways of computing f,(x) are given by the following methods 


I, I] and IIT. 





Method I. The obvious way to compute f,(x), if 7 is given, is to carry 
out the successive additions and divisions indicated by the form of the 


TaB_eE 1. Computation of arc tan 1 using ten terms of the continued fraction 


i bys , haut a 

0 19 19. 81 4.26315789 
1 17 21.26315789 64 3.00990099 
2 15 18.00990099 49 2.72072567 
3 13 15.72072567 36 2.28997063 
4 11 13.28997063 25 1.88111778 
5 9 10.88111778 16 1.47043717 
6 7 8.47043717 9 1.06251895 
7 5 6.06251895 4 .65979175 
8 3 3.65979175 1 .27323959 
9 1 1.27323959 1 -78539814 


expansion when it is written as in (1). More formally, the method consists 
of calculating the sequence 


—= Bans + Cn i+1 Cn+1 = 0 
Co = An-i/dn 


fori = 0,1, --- m — 1. Then f,(x) = bo + c. This method is illustrated in 
Table 1 by the computation of arc tan 1 using 10 terms of expansion (5). 
Method II. Successive approximants, for consecutive n, may be calcu- 
lated by writing 
fir(x) = An/ Ba, 


where 


Ay 
Bo 


bo A, = bob; + a; Ant = bnyiAn + OnytAn—s 
1 B, = b, Basi = bn Bn + Anyi Bri. 


An example using this method is given in Table 2 where expansion (4) is 
used to calculate In 2.3026. 
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TABLE 2. Computation of In 2.3026 by Method II 





n bn a, Ay B, A,/B, 
0 0 0 sl xX 10! 
1 1 1.3026 .13026 x 10! 1 xX 10! 1.3026 
2 2 1.3026 .26052 X 10! .33026 xX 10! -788833041 
3 3 1.3026 .951236676 X 10! .112104 xX 10 -848530539 
4 4 5.2104 .5162360112 K 10° .62049467 xX 10° .831974933 
5 5 5.2104 .3076812414 x 10° .3686580032 * 10° -834598025 
6 6 11.7234 .2451291574 x 10' .2939378741 X 10# .833948868 
7 7 11.7234 .2076611128 « 10° .2489757642 « 10° -834061554 
8 8 20.8416 .2172177286 X 10° .2604419673 x 10° .834035048 
9 9 20.8416 .2387758543 X 10° .2862883034 X 10° .834039852 
10 10 32.5650 .3095128077 « 10° .3711012300 x 108 -834038754 
11 11 32.5650 4182214453 x 10° .5014411389 10° .834038959 
12 12 46.8936 6470074325 x 10 —_.7757520931 x 10” 834038913 
tai 13 13 46.8936 1037228754 X 10" 1243621523 X 10" 834038922 
; 14 14 63.8274 1865088276 X 10" — .2236212522 x 10" (834038919 
15 15 63.8274 3459668557 X 10" 4148090065 X 10" ‘834038920 
Method III. The nth approximant may also be expressed as a sum, i.e., 
hods frlx) = bo + D pipe --> pi, 
i=l 
where 
wd i ite 1 
2 Tr; = , —a i = a 
- bib; 1+ rd — wd 
with initial values 
ay, 1 ‘i 1 
— b vail 1+r.° 
Table 3 gives the computation of In 2.3026 by this method. 
TABLE 3. Computation of In 2.3026 by Method III 
n Tn 1 + pp Pip2 *** Pn f 
1 1.3026 1.0 1.3026 1.3026 
2 6513 .605583480 — .513766959 -788833041 
3 2171 .883804325 .059697499 -848530540 
4 4342 .722675023 — .016555608 .831974932 
5 .26052 -841558660 .002623093 -834598025 
sists 6 .39078 .752522289 — .000649157 -833948868 
7 .279128571 .826411810 .000112686 .834061554 
8 .372171429 .764779246 — .000026506 -834035048 
9 .289466667 .818747283 .000004804 -834039852 
10 361833333 .771456089 — .000001098 .834038754 
11 .296045455 .814076311 .000000204 .834038958 
12 .355254545 -775672245 — .000000046 .834038912 
d in 13 .300600000 .810920125 .000000009 -834038921 
(5) 14 .350700000 .778579904 — .000000002 -834038919 
| 3 -“ .303940000 .808641743 .000000000 .834038919 
Icu- 


3. Truncation Error. These errors arise because only a finite number 
of terms of an infinité process can be used. For a discussion of these errors 
it is convenient to separate functions into classes on the basis of the 
order of magnitude of the mth term in their power series expansions. If the 
order of magnitude of the mth term is x"/n! the function belongs to one 
class; if it isx"/n the function belongs to a second class, etc. The reason for 
1) is this classification is that functions in the first class can generally be com- 
puted adequately from their power series, while for functions of the second 
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class an unreasonable number of terms may be required. HEISING” states 
that in computing arc tan x by means of a general purpose board on the 
604 Electronic Calculating Punch the maximum error after 990 terms are 
used is still 11 X 10-7. Fortunately the continued fractions of some of the 
functions in the second class converge rapidly enough to make them prac- 
tical for machine computation. 


4) 
poe) eo 


7 

. . f 

TABLE 4. Truncation error in e* after 10 terms re 

x Power Series Continued Fraction t 

1. 2.7 X 10% .67 X 10°% ] 
2 61 X 10-* 22 X 10°* 
3 59 X 1074 44 X 10-4 


Exampies of functions belonging to the first class are e’, sin x, cos x, 
sinh x and cosh x and examples of functions belonging to the second class 
are In x, arc sin x, arc sinh x, arc tan x and arc tanh x. To compare the con- 
vergence of the power series and the continued fractions expansions of | 
functions in these two classes we have chosen e” to represent the first class 
and In x and arc tan x to represent the second. The continued fractions 
expansions are given by (3), (4) and (5) and the power series are | 





. 2 3 
(7) e=1+5+5+5t- 
— 1 — 1)? Lae am. Op 
(8) inx = &= :.2 > c42 3 er 


9) ; - 2. 

arc tanx = = 
‘ a ae 
Tables 4, 5 and 6 give the comparisons. The calculations were carried out 
to two more decimal digits than were required, this should make the round- 
ing error negligible (except, perhaps, for the computation of In x for x near 


zero). 


TABLE 5. Number of terms required to compute In x to 9 decimals 


x Power Series Continued Fraction 
.0001 95,000* 550 
.0010 11,500* 315 
.0101 1,350* 105 
.1054 143 34 
5108 25 16 
.6931 16 11 
.9163 8 7 

2.3026 — 16 
4.6052 — 24 
6.9078 _ 30 


* Approximate values. 


Table 4 gives the difference between e* and the 10th approximant from 
the power series and from the continued fraction (in the power series the 
last term is x/10!). This table indicates that the continued fraction con- 
verges a little faster than the power series but the difference is too small to 
be of practical importance. 
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Table 5 gives a comparison of the number of terms required for fixed 
accuracy. For the power series the number of terms, 1, is determined by the 
condition that the absolute value of the (m + 1)th term in (8) shall be less 
than 9.3132 X 10-". The continued fractions were computed by Method III 
and n was obtained as the first integer for which pipe --- px < 9.3132 & 10-”. 
The table shows that the continued fraction expansion converges appreciably 
faster than the power series in 0 < x < 2. Actually the continued fraction 
converges for all x > 0, while the power series diverges for x > 2. The con- 
tinued fraction could be used to compute, in a reasonable number of terms, 
In x for .1 <x < 10. 

Table 6 gives the corresponding comparison for the arc tangent. For the 
power series, » is determined by the condition that the absolute value of 
the (n + 1)th term in (9) shall be less than 10-* while for the continued 
fraction n is the number of terms required to reduce the truncation error 
to less than 10-*. The continued fraction converges faster than the power 
series and also converges in the region where the power series does not, 
i.e., for x greater than 1. 


TABLE 6. Number of terms required to compute arc tan x to 6 decimals 


x Power Series Continued Fraction 
1 3 3 
Pe 4 3 
e 5 4 
A 6 5 
Ss 8 5 
6 11 7 
a 15 7 
38 22 7 
9 44 8 
1.0 - 8 
2.0 15 


4. Machine Computation. Computing machines differ greatly in many 
aspects. In deciding what method should be used to compute values of 
functions the amount of high speed storage available is one of the most 
important factors. ’ 

Machines which have enough high speed storage to accommodate arbi- 
trary coefficients can probably compute values of many functions most 
rapidly by using polynomial approximations with a fixed number of terms. 
If, however, the majority of the arguments fall in the region in which only 
a few terms of an infinite process are required, or if varying degrees of 
accuracy are required from one computation to another, the infinite process 
may be more efficient. In general the computation of a term in a power 
series requires fewer operations than the computation of a term in the 
continued fraction. However for functions such as In x and arc tan x more 
terms of the power series are required. In either case it is desirable to com- 
pute only as many terms as are required to achieve the required accuracy. 
Therefore, for computing continued fractions, Method I would not be 
suitable, but Method II or III could be used. It is immediately evident 
from Table 2 that Method II cannot, in general, be used without a floating 
decimal point, because A, and B, increase too rapidly. Furthermore, since 
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Method II requires as many operations and more storage than Method III 
the latter is more suitable for machine calculations. 

Table 7 gives the calculations required in computing successive approxi- 
mants to the continued fraction, the values required from storage and the 
initial values, if Method III is used. If a, and 6, are simple functions of n 
and x which can be computed as required, the first calculation requires the 
storage of only m and x. For example, for expansion (4), 


i, = [3n P(x mr 1)/{n(n a 1)}, 


where [yy] denotes the greatest integer <y. Even with this simplification, 
five quantities must be stored. 

This fact makes this method of computation impossible for machines 
with very small internal memories, such as the IBM 604 Electronic Calcu- 
lating Punch. The computation of power series requires very little storage 
and can be carried out on such machines. However for functions of the 
second class where too many terms of the power series are required, con- 
tinued fractions can be used if the maximum number of terms of the 
continued fraction expansion that will be needed to give the desired accuracy 
is determined in advance and Method | is used to compute the value of the 
expansion. If the coefficients can be generated in a simple manner only three 
numbers have to be stored at any one time."! 


TABLE 7 
Calculation Required from Storage Initial Values 
&. fn = An/dnrbn — =< % 1 
2. 1 + Pn = {1 +r,(1 + Pn—1)}~* 1 + pn-1 1 
3. (pip2 *** Pn—1)Pn Pip2 *** Pn a/b; 
4, Sn = fn t+ pip2 *** pn Ju~a bo + ai/hi 


5. Bibliography. The general theory of continued fractions is given by 
PERRON! and WALL.’ Both of these books have extensive bibliographies. 
The use of continued fractions in interpolation is discussed in MILNE- 
THompson,’® NORLUND‘ and LANE.® Mutr® gives some general methods for 
transforming infinite series into continued fractions. MULLER,’ AROIAN® and 
BurGEss° give some examples of the use of continued fractions in evaluating 
integrals. 

D. TEICHROEW 
University of North Carolina 
Raleigh, N. C. 
NBSINA 


This work was performed on a National Bureau of Standards contract with the Uni- 
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5R. E. Lane, “Interpolation by means of continued fractions,” Fraternal Actuarial 
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7 J. H. MOLLER, “On the application of continued fractions to the evaluation of certain 
integrals, es special reference to the incomplete Beta function,” Biometrika, v. 22, 1920-1, 
p. 284-297. 

8 L. A. Aroran, “Continued fractions for the incomplete Beta function,” Annals Math. 
Stat., v. 12, 1941, p. 218-23. 


® J. Burcess, “On the definite integral =f" e~® dt with extended tables of values,” 


Roy. Soc. of Edin., Trans., v. 39, part II, 1898, p. 257-321. 

10W. P. Hetsinc, “An eight-digit general purpose control panel,’”” IBM Technical 
Newsletter, no. 3, 1951. 

1 The use of continued fractions in computing In x, arc tan x and arc sin x on the 
Model II Card Programmed Calculator will be discussed in a forthcoming National Bureau 
of Standards report. 


On a Punched-Card Method of Solving 
Certain Integral Equations 


1. Introduction. In the present paper we are concerned with the numeri- 
cal solution of the homogeneous Fredholm equation 


1 
(1.1) 660) = [ ote, vote) dx, 
where the given kernel p(x, y) ‘satisfies the conditions 
(1.2) p(x,y)>0, OSx<1, OSy<1 
and 
1 
(1.3) [onda t, oss. 
0 


The solution ¢(y) is to be a non-negative function such that 


(1.4) [eora = 1. 


Equations of the type (1.1) have received considerable attention in the 
theory of probability,!? and the above problem has been formulated with 
this in mind. Although all the examples which have been considered by the 
authors are special cases of (1.1), the numerical method to be explained is 
immediately applicable to the following more general equation 


b 
(1.5) $(y) = fly) + f K(x, y, \)6(x) dx, 


where f(y), K(x, y,\) are given functions, \ is a parameter (supposed 
known) and a, b are constants. This integral equation includes the general 
non-homogeneous Fredholm equation 


(1.6) 60) = foo) +r f° K(x, 9)o(@) de 


as a special case. 

Numerical methods for the solution of integral equations have been 
considered by a number of writers." The methods employed may be 
classified briefly as follows: (a) Methods which involve the solution of a 
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system of simultaneous linear equations, (b) Iterative procedures, (c) Re- 
laxation methods, (d) Monte Carlo methods. Of these, only (b) has been 
found suitable for use with punched cards and such a method of solving 
equation (1.1) is explained below. 

2. The Method of Solution. It is clear from (1.4) that the function ¢(y) 
may be regarded as a probability density function. It is then natural to 
approximate its graph by means of a histogram. The general form of this 
histogram will often be suggested by the particular problem or may be 
ascertained by a study of the kernel p(x, y). Thus, as a first approximation, 
we replace ¢(y) by the step-function 


(2.1) dy) = 6, (f-D/nSy<j/m, GF =1,2,-+0-1, 
=c,, (n—1)/n<y <1, 


where 7 is a positive integer—the number of steps in the histogram—and 
the c’s are a set of positive constants. It is not essential at this point that 
these constants should satisfy (1.4), i.e., that 


(2.2) pa cj = Nn. 
j=1 
Substituting from (2.1) in the right member of (1.1), we find an approxi- 
mation to ¢(y), viz. ¢*(y), where 


n 


a 
(2.3) ¢*(y) =X a0 p(x,y)dx, OSy<1. 
(j-1)/n 


j=l 


Since p(x, y) is known, the integral appearing in the above equation may be 
evaluated for any given value of y. Let the value of this integral be denoted 
by A,(y). Then equation (2.3) becomes 


(2.4) ¢*(y) = LD c/Aj(y), O<y<l. 
7=1 


In particular, (2.4) provides an approximation to ¢(y) at each value y = 
(k = 1,2,---,m), where », represents an arbitrary point in the kth inter- 
val ((k — 1)/n < y < k/n). In practice, it was found convenient to take 
ye = (k — 1)/n. The values $*(yx) (R = 1, 2, ---, 2), calculated from (2.4), 
provide a new set of constants c,“ which in turn can be used to define a 
new step function ¢;(y) to replace ¢o(y). Explicitly, 


(2.5) ch? = Do cjAj(y«), k= 1,2, ---,m, 
jal 
and 
_ ja®, (k—1)/n<y<k/n 
(2.6) oi(y) = - (n —1)/n<y <1. 


If necessary, the constants c,“” may be normalized to satisfy (1.4) by 
multiplying each c. by n/( ; <j”) . A further approximation ¢2(y), de- 
j=1 


fined in terms of constants c,, can be found in precisely the same way, and 
this procedure may be continued indefinitely. 
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This iterative procedure is continued until the normalized values of the 
constants ¢,“ from two successive iterations, r = m and r = m + 1, differ 
by less than some pre-assigned amount. An analytical proof that this pro- 
cedure converges is given by BERNIER.‘ In practice, the final approximation 
depends somewhat on the number n of steps in the approximating histo- 
gram—the larger m is, the better the approximation. 

3. The Plan of Computation for Punched Cards. The numerical calcu- 
lations involved in the iterative procedure described above may be per- 
formed using punched card equipment. The following IBM machines, or 
their equivalent, are required : 


1. The Alphabetical Accounting Machine, Type 416, 
2. The Calculating Punch, Type 602-A, 

3. The Reproducing Punch, Type 513, 

4. The Sorter, Type 070. 


An outline of the plan for carrying out these calculations is given below. 


I. Prepare a deck of n* cards, to be called the A deck. On the [(k—1)n+j ]th 
card, (k, 7 = 1, 2, ---,), punch the following information: 


(i) the card number (k — 1)m + 7 in card columns 78-80, 
(ii) the numerical value of A ;(y,) in card columns.1-5, 
(iii) the group number j in card columns 75, 76. 


The number of card columns required for each of the fields (i) to (iii) will 
depend on the particular problem; the positions of these fields on the card 
are, of course, arbitrary. 

II. Prepare a second deck of cards, to be called the C deck. On the jth 
card (j = 1, 2, ---,), punch the following information: 


(i) the card number j in card columns 75, 76, 
(ii) the number jn in card columns 78-80, 
(iii) the number c; in card columns 65-70, 
(iv) an X in card column 7. 


For succeeding iterations the numbers c; in field (iii) will be replaced by 
cs (r = 1,2, ---). Field (iv) serves a dual purpose; it provides a ready 
means of distinguishing the cards of the C deck from those of the A deck 
and it also provides the means of controlling multiplication (cf. III below). 


III. Carry out the following steps in order: 


(1) Place the C deck in front of the A deck and sort on card columns 
75, 76. The jth card of the C deck, bearing the number c;, is now followed 
by all those A cards in which are punched the numbers A (yx) to be multi- 
plied by c;. 

(2) Using group multiplication, form the products c;A (yx) and punch 
each answer in card columns 8-15 of the corresponding A card. 

(3) Sort out the C cards, place them behind the A deck and sort on 
card columns 78-80. The A cards are now arranged in m groups, the kth 
group being followed by the kth C card. The sum of the products ¢;A ;(yx) 
in each group will give the corresponding unnormalized value of cj. The 
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cards of the C deck are here used merely to provide a break in the control 
and so allow a total to be taken. 

(4) Add the products c;A,(y,) and list the resulting sums c. This 
step is performed on the tabulator and the next deck of C cards containing 
cx in place of c, may be prepared simultaneously by summary punching. 
(5) Sort out the old C deck. 


The first iteration is now completed. 
It is not necessary to calculate normalized constants c,“? at each stage, 
i.e., a set of constants satisfying 


(3.1) La” =1 

k=1 
However, since the criterion for terminating the procedure depends on a 
comparison of normalized constants, we must eventually find such nor- 
malized values. When they are required, these are easily found. 

4. Illustrative Example. Consider a point moving at random on a 
segment of the x-axis lying between reflecting barriers at x = 0 and x = 1. 
Let p(x, y) dx denote the probability density of the point moving in one 
step from a position in the interval (x, x + dx) to a position with abscissa y. 
Then ¢(y) dy is the probability of the point being in the interval (y, y + dy) 
after an infinite number of steps. 

With regard to the actual motion of the point, we make the following 
assumptions : 


(a) The point moves in steps whose lengths are normally distributed 
with mean .1 and standard deviation .02. It follows that single steps involv- 
ing more than one reflection at a barrier have such small probabilities that 
they may be ignored. 

(b) Steps from the central reference point yo = .5 take place in either 
direction with equal probabilities. 

(c) Starting from any point other than yo, steps toward yo have proba- 
bility .7 whereas steps away from yo have probability .3. 


It follows that 





.7g(y—x) +03 gyt+x), OSxeSy¥<y, 
3g(x—y +03 g(e+y), OS¥ySx<y, 
.7 g(x — y), y¥<% xX > Yo Or y> Yo 

(4.1) p(x,y) = + x<yoOry=y x HF Yo, 
3g(y—x) +03¢(2-x-y), wexSySl, 
Tg(x—y) +03 g(2-x-y), weysxSl, 
[Sg(x—y), «=H OSy<1, 

where 

(4.2) g(t) = (2x)-4o— exp { —3(t — .1)?/0?}  (o¢ = .02). 


From a consideration of the assumptions (a) to (c) or from a study of a 
graph of p(x, y), we easily see that ¢(y) must be symmetrical about the 
point y = 0.5. Thus, if the number of histogram steps is even, we must 
have 


(4.3) co” = OP seas y = 0,1,2,---. 
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Fic. 1. Limiting n-step histograms approximating ¢(y), » = 10, 20, 40, 80. 


Here, as elsewhere, r represents the number of the iteration. This means 
that we need only compute the ordinates of ¢;-)(y) over one-half the range 
(0, 1). 

Approximate results for lim ¢)(y) were obtained for m = 10, 20, 40 


and 80. In order to obtain agreement between successive iterations to three 
significant figures, the numbers of iterations required were 12, 24, 33 and 32, 
respectively. Figure 1 shows a graph of the four limiting histograms, and 
it is clear that there is quite good agreement among them—they appear to 
be approximations to the same function. The numerical convergence is 
exemplified by the successive normalized values of c; shown in Table I for 
the case m = 10. 


TABLE 1. Normalized values of c; for n = 10 
r = the number of the iteration; c;" is the value of $)(y) 
in the interval ((j — 1)/m < y < j/n). 


x 1 2 3 4 5 6 
j 

5, 6 2.18 1.84 2.05 2.05 2.01 2.00 
4, 7 1.57 1.33 , 1.50 1.24 1.50 1.47 
3, 8 0.73 0.67 0.80 0.74 0.80 0.79 
2, 9 0.38 0.44 0.45 0.47 0.47 0.49 
1,10 0.19 0.07 0.19 0.23 0.26 0.26 
\r 7 8 9 10 11 12 
j 

5, 6 1.93 1.96 1.97 1.97 1.97 1.97 
4, 7 1.48 1.45 1.46 1.46 1.46 1.46 
3, 8 0.80 0.80 0.80 0.80 0.80 0.80 
2, 9 0.51 0.51 0.51 0.50 0.50 0.50 
1, 10 0.28 0.28 0.28 0.28 0.27 0.27 
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In this example, the time required per iteration was‘about 15 minutes 
for the case m = 20. For the larger values of 1, the time required was only 
slightly longer. Since convergence was fairly rapid, a good approximation 
to the solution ¢(y) could be obtained in a few hours, once the initial A 
and C decks have been punched. 

The method has been applied to about eighteen different cases of the 
homogeneous equation (1.1). The solutions obtained agree well with those 
found by Monte Carlo methods and are supported by a limited body of 
experimental evidence resulting from a study of the one-dimensional move- 
ments of fish. 

B. A. GRIFFITH 


K. W. SMILLIE 
Computation Centre 
University of Toronto 
Toronto, Ontario 


Part of the work described above was sponsored by the Ontario Research Council. 
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The Use of Exponential Sums in Step by 
Step Integration—II 


[Continued from MTAC, v. 6, p. 63-78] 


10. The error analysis for the exponential method is predicated on the 
convergence of the series representation of the function 


Wiitynait+y}, yoe*-1 
about y = 0 (see equation (13)). The radius of convergence of this series 
is 1. For small h, y = e~” — 1 is small and there will be convergence for all 
positive h less than the least h = fo such that 
(32) jem — 1] = 


Such an ho exists except when A is real and non-negative, in which case 
|y| < 1 for every positive h. 
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For rapid convergence and efficient use of the method, it is imperative 
that 4 be much smaller than ho. If X = a + 48 is not real and non-negative, 
there is a simple way of computing fo based on a table obtained from 
equation (32). 

Equation (32) is equivalent to 


(33) e-* = 2 cos Bhp. 


Let u = aho, v = Bho. For a ¥ 0, let B/a = v/u = R. If a = 0, Bho = 2/3. 
This is of interest in itself. We may choose » to be in the first quadrant. 
(34) u = — In (2 cos?). 


For values of v between 0 and 2/2, u may be determined, and hence, the 
ratio R. To use the table, determine R as the ratio B/a, B > 0. The table 
then gives the corresponding u and v from which hy may be determined 
either as 


(35) ho = u/a or ho = 0/8. 


Table of u and R as functions of 9 


v u = —In (2 cos») R=v/u 
0 — .69315 0.00000 
R — .68814 —.14532 
a — .67301 —.29717 
a — .64746 — 46335 
4 — .61092 — .65475 
ss — .56256 — .88879 
6 — .50118 — 1.19717 
7 — 42506 — 1.64682 
8 —.33176 — 2.41141 
9 —.21770 — 4.13404 
1.00 —.07752 — 12.89978 
1.01 — .06177 — 16.35005 
1.02 — .04567 —22.33275 
1.03 — .02921 — 35.26553 
1.04 — .01236 — 84.11687 
1.05 .00487 215.61735 
1.06 .02251 47.09604 
1.07 .04056 26.37861 
1.08 .05905 18.28862 
1.09 .07799 13.97557 
1.1 .09740 11.29353 
1.2 .32198 3.72698 
1.3 .62549 2.07836 
1.4 1.07900 1.29749 
1.5 1.95564 -76701 
1.51 2.10770 -71642 
1.52 2.28721 66456 
1.53 2.50629 -61046 
1.54 2.78737 55249 
1.55 3.17990 48744 
1.56 3.83542 40673 
1.57 (/ 644235 .24370 


11. The calculation of the coefficients, ao, «++, @n-1, of the open method 
of integration may run into considerable difficulty because of the small 
determinant involved. In this section practical methods for computing 
these coefficients without great losses of accuracy are given. 

The coefficients, ao, @:, -**, @n-1, for the open method of integration 
are determined by equations (4). By means of the quantities x; = «~”", 
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3 = x; — 1 and the function 





Yi 
36 ) = 
(36) £0) = Ty) Ind +) 
described in § 8, (4) may be written 


n—l 


(37) 2 a,x = g(yj). 


One can then evaluate the a;’s by Cramer’s rule. Thus, if @ is the Van- 
dermonde determinant of the x;’s and a; the result of replacing the ¢ + 1 
column by g(y;), 


(38) at = a,/a. 


This can be evaluated directly. However, to obtain control over the accu- 
racy required for these coefficients, it is desirable to arrange the computation 
so that the loss in accuracy in computing these coefficients is a minimum. 
The following is believed to be a practical method for this purpose. 

Let 


n—1 


(39) g(yi) = Do Any* + hs. 
k=0 


We will manipulate the determinant a, by its columns only and we will 
indicate it by a row without subscripts, i.e., 


(40) ar= 1, ees, g(y), en, ooo, get], 
As in the discussion of «(A/%) following equation (24), one may replace 
x; by y; in the first ¢ columns. 


n—1 


a, = |1,---, yy, ¥ Any? + AM, xt, ---, x], 
k=0 


Now x = 1 + y and by an obvious reduction one obtains 
n—-! 


1, si “7 > Axy* + hk”, yy" + (t + 1)y'‘, sedan rs 


k=t 





ayg= 


- 1 
+ (n — lyr? +. +("7 Jor 





By an elementary argument, the latter columns in this determinant can 
be reduced so that 


n—l 
> 1. a “+. an, } Ary" + h”, : ae + (é + 1)y'‘, it yy 
k=t 


+ceon('t 4p, 


= r oe ~~, Dy' ote h™, y* + (é + 1)y'*, vey yit* 


+(- (Fy, |, 


in- 
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ion 
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where 


t+2 
(41) Dim Ac~ (+ Am + ("9 Aue + 


+(- yr(” P Aes 
Thus, 
a: = Dat It see, yt, AM, yt + (t + I)y', ---, yt 


+ (= y(’ : ys, vl, 
Let 


t+k 
B: = t soe, yh hm, yt + (¢ + 1)y', ren, ytth + (— 1y( : Jy | 
= L (- 1)*+#+14 6 | 1, ae + y*', yt, = °, 77" | ne j 200 
7 
+ (t + 1)|1, coon gr, y', y, +20, 9" | n0 § row 
t+2 
+( r Jit, 9 a ot, 22+, "| no j row + ssiadlins 


nm—i1- 
bs . —t-— Jit “+5 9" | no § row]. 
Now, using the notation of § 7, we have 
; n—t—1 t oe k 
B: = LX (= 1yererno( } * ( t We tek (Vay °° *» Vinay Visas °° +9) 
7 k=0 
n—t—1 t 
= om (-— 1yeeetngo( > ( **\- 1)***V, (no yi) 
i k=0 
x (- 1) ***D¢n-1)—-+-2S0, +) 


a2 = E(—no(S"('F*)v0 (no »bee-n-r4 (n0 9.) 


i k=0 t 


= x (— 1)*42,™ Vy (no w(E (" ‘ Vous (no »)) 


k=1 
—! /(n —k 
= a(— 1) HD Ai (resi — yA E ( , ons (no »)) 
? 
= a(— 1)*H(>> KE; t). 
i 
Thus, one can calculate a, by the formula 
a: = a/a = De + (— 1)*4 KE; 
i 
where D, is given by (41) 
(43) Ki = hi /TWies(¥5 — 9) 


(44) Exe= 2 (” ’ "Vs (no yj) = = e cs . rhe pn (no 4). 


r=x0 











142 USE OF EXPONENTIAL SUMS IN STEP BY STEP INTEGRATION 


It is of interest to note that D, is the coefficient corresponding to the 
case in which each »; is zero. 

The critical stage in the computation of a; is in computing h;™ for K;™, 
because of the cancellation of many significant digits. 


n—l 


hs = g(yi) — LD Any#. 
k=0 
It is probably best to use u; = hv; (h is the time interval) as the basis 
for this computation, since the »; are chosen precisely. Thus, 
y= evi— 1, gyi) = (C1 — 1)/u;. 
The A;’s are most readily computed by the inductive formula 


Oe ties) dette Malt 2 @1,2,---. 


This formula is obtained from the fact that g(y) = Ao + Ary + Asy? + --- 
is the reciprocal of (y + 1) log (1 + y)/y whose series can be readily written 
out. The loss of accuracy can be judged from the computation of fh, in 
§ 12 below. Thus, for nine places after the decimal point: 


g(v:) = .98457 9722 + .11180 85987 
5 
X Aun* = .98462 3553 + .11180 4513 7 
k=0 

h,® = —,00004 3831 + .00000 4081 i. 


The following computational procedure is suggested for computing the a,. 


1) The »;, h and m are chosen and the u; = hy; calculated. 

2) The D; are obtained. 

3) From the u;, y; = e~“i — 1 and g(y,;) = (e*i — 1)/u; are obtained 
and used in turn to compute /;™. 

4) The various differences y; — y; are computed and from these the 
products Ili«;(y; — y;) are formed. These can be used to compute K;™ in 
conjunction with the result of 3. 

5) The ~; (no y;) are obtained by forming the polynomial 


y=! dia py" + poy" 4. + oe (- 1)*"Do. 


The complex y; are combined in conjugate pairs to yield the quadratic 
factors for this polynomial except, of course, for the conjugate of the y; 
omitted. The ~; are then used to obtain the Ej». 

6) a; is then obtained from the above formula. 

The process 3) is readily seen to be the only one critical as far as accu- 
racy is concerned and overall estimates for the accuracy of a; can be obtained 
based on this. 

Let a’, ---,a8_; be the result of the above computation. When one 
substitutes back in (4) or (37) one may find a discrepancy e; in each equa- 
tion, i.e., 

n—1 


a a,x 7 = gyi) — €;. 
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Using (37) one obtains 
n—1 
xX @ — a))xf = €j. 
raf) 


Thus if a, = a,° + Aa,, Aa, satisfies an equation 


n—l 


> Aa,x7 = €;. 
r= 


Now the previous argument can be used to show that 
(46) Aa, = (— 1)***1 5 LE; t, 
i 


where L; = ¢;/Ili«;(y; — ys) and the E;,, are the quantities previously com- 
puted. Thus, the computation can be readily continued as far as the accuracy 
to which y; and g(y;) will justify. 

12. A report on the aspects of numerical solutions of differential equa- 
tions pertinent to analogue computing installations is being prepared by 
Project Cyclone. As part of the report, a system of differential equations 
representative of the types encountered is solved by various numerical 
procedures. The exponential method gave the most accurate solution, 
although it was used with a very large interval width. The details of the 
problem as solved by the exponential method are considered in this section. 

The problem is the following: 


V = 9.295 cos a — 32.2 sin y — .00056022 V°[.129 + .051632(.965 + 5.1a)?] 


V¥ = 9.295 sin a — 32.2 cos y + .00056022 V*(.965 + 5.1a) 
y = —.00009421 V*(.215y + .44a — .026) 


6=y 
a=0-—y¥y. 
The initial conditions of the system are: 
V(0) = 200 y(0) = —.0204 
7(0) = 0 6(0) = .0525 


Solutions of the system of equations had been obtained by other pro- 
cedures before the exponential method was used. Based upon these solu- 
tions, characteristic matrices were calculated at ¢ = 0 and 5.7 seconds. The 
characteristic values were 


t AL re As, M 

0 —.72140 2212 + 1.28266 5342 — .01586 3538 + .19708 34222 
5.7 — .86065 6195 + 1.43671 6492 — .01938 6195 + .17749 50264 
In the exponential method m = 4 and », v2 = —.80 + 11.36 and vs, % 


= —.018 + .19i. These values are rounded averages obtained from ),, As 
and As, As, respectively, of the above table. It is clear from the table of 
values that the chosen quantities are not very different from the charac- 
teristic values of the ¢ = 0 matrix that would normally be chosen. 

Based on the characteristic values, the short period of the solution is 
approximately 4 seconds. An # value was chosen by a rule generally applied 
to cubic polynomial approximation methods, namely, h is not to exceed 
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T/12. Thus, h was originally chosen as .3. The natural limitation on the 
size of h, discussed in § 10, shows that h must be chosen less than .9. This 
was obtained by evaluating a maximum h for the characteristic values of 
the ¢ = 0 and ¢ = 5.7 second matrices. 


t a B max h 
0 —.721 1.28 1.041 
5.7 — .861 1.44 .934 
0 — .016 .197 5.562 
52 —.019 177 6.276 


Having chosen an interval width, it is necessary to estimate the error of 
the solution. A-bound on the error can be obtained by an evaluation of the 
step error formula, equation (10), using the relationships obtained for 
e(A;h) and e(t, c) in § 8, 9. An estimate of the total error in the solution was 
obtained by a recomputation at a half interval width, h = .15. Since the 

' step error is proportional to h', the truncation error of the second procedure 


TABLE A 
t V y 
h= 3 hk = .15 A108 h= 3 h = 15 410° 
0 200.000 200.000 .00000000 .00000000 
200.700 — .00298124 
3 201.407 201.407 — .00576910 — .00576886 
202.121 — .00833298 
6 202.840 202.840 — .0106580 — .0106575 
203.562 — .0127295 
9 204.286 204.286 — 0145486 — .0145079 
205.010 — .0161168 
1.2 205.734 205.734 0 — .0174438 — .0174457 19 
206.457 — .0185476 
1.5 207.178 207.178 0 — .0194343 — .0194384 41 
207.897 — .0201350 
1.8 208.612 208.612 0 — .0206523 — .0206549 26 
209.325 — .0210150 
2.1 210.033 210.033 0 — .0212278 — .0212311 33 
210.738 — .0213174 
2.4 211.438 211.438 0 — .0212849 — .0212861 12 
212.134 — .0211476 
2.7 212.824 212.824 0 — .0209110 — .0209099 —i1 
213.510 — .0205790 
3.0 214.190 214.190 0 — .0201594 — .0201593 -1 
214.864 — .0196533 
3.3 215.532 215.532 0 — .0190644 — .0190624 —20 
216.193 — .0183866 
3.6 216.847 216.847 0 — .0176267 — .0176256 —11 
217.494 — .0167781 
3.9 218.132 218.132 0 — .0158419 — .0158429 10 
218.762 — .0148186 
4.2 219.383 219.383 0 — .0137044 — .0137039 - 5 
219.995 — .0124979 
4.5 220.596 220.597 —1 — .0111973 — .0112000 17 
221.188 — .00980993 
4.8 221.768 221.769 —1 — .00832525 — .00832802 26 
222.338 — .00675493 
5.1 222.896 222.896 0 — .00509068 — .00509185 12 
223.442 — .00334027 
5.4 223.975 223.975 0 — 00149742 — .00150204 46 
224.496 .000420690 
5.7 225.004 225.004 0 .00242675 .00242555 12 
225.499 .00451008 
6.0 225.980 225.980 0 .00667289 .00667159 13 
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TaBLe B 
t y 4 

k= 3 h = .15 410° h=.3 h = 15 410° 
0 — .0204000 — .0204000 .0525000 .0525000 
— .0165274 .0497366 
B — .0131725 — .0131726 0475157 .0475158 
—.0103520 .0457580 
6 — .00805566 — .00805626 0443837 .0443841 
— .00625330 .0433168 
9 — .00489427 — .00489545 .0424856  .0424860 
— .00392390 .0418290 

1.2 — .00328144 — .00327280 — 86 .0412954  .0412928 26 
— .00287426 .0408345 

1.5 — .00267128 — .00266115 —101 .0404182 .0404212 — 30 
— .00257054 .0400300 

1.8 — .00254925 — .00254868 — 6 .0396415  .0396467 — 52 
— .00253762 .0392653 

2.1 — .00250967 — .00250611 — 36 .0388805 .0388865 — 60 
— .00242001 .0385161 

2.4 — .00225020 — .00225703 68 .0381541 .0381643 -—102 
— .00200326 .0378435 

2.7 — .00164177 — .00165228 105 .0375630 .0375681 — Si 
— .00120409 .0373527 

3.0 — .000656351 —.000664002 76 .0372078 .0372115 — 37 
— .0000412040 .0371577 

3.3 .000665564 .000652188 134 .0372004 0372027 — 23 
.00140242 .0373562 

3.6 .00219999 .00219509 49 .0376306 .0376256 50 
.00301588 .0380042 

3.9 .00385298 .00385141 16 .0385338  .0385311 27 
.00468956 .0391718 

4.2 .00552401 .00552000 40 .0399427 .0399376 51 
.00633420 .0408270 

4.5 .00711795 .00712560 — 76 .0418444 .0418368 76 
.00788956 .0429633 

4.8 .00862167 .00862310 — 14 .0442031 .0442021 10 
: 00932485 .0455486 

5.1 .00999301 .00999468 -— 17 .0470037 .0469980 57 
.0106334 .0485455 

5.4 .0112334 .0112428 — 94 .0501893 .0501865 28 
.0118249 .0519169 

5.7 .0123873 .0123822 51 .0537306 .0537328 — 22 
.0129170 .0556305 

6.0 .0134275 .0134319 — 44 .0576129  .0576069 60 


may be considered negligible in comparison with that of the first compu- 
tation. The difference in values of the solution at h = .3 and h = .15 was 
very small. These values are given in Table A. The maximum errors given 
as a percentage of the absolute maximum of the variables V and y are .000 
and .022, respectively. 

By the use of § 10 to calculate the coefficients of the integration for- 
mulas, the following values were obtained: 


h ao a ae a3 
a 2.09056 049 — 1.92175 461 1.07630 795 — .23931 8630 
15 2.19751 098 — 2.19425 286 1.29578 949 — .29906 1942 


The corresponding results for the variables y and @ are given in Table B. 
The percentage errors in these cases are .066 for y and .020 for @. 

The initial four values to start the procedure were calculated in the 
h = .15 case using the Runge-Kutta method. In the h = .3 case these values 
were calculated using an unpublished method of R. F. CLIPPINGER. 
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13. The interest of the authors in the use of exponential sums for inte- 
grating differential equations arose from their experience in integrating a 
system of fourteen first order differential equations in order to provide a 
large scale test problem for a continuous computer. For references to this 
work see MTAC, v. 6, p. 78, ftn. 3. 

The system is stable and for each of the fourteen unknown dependent 
variables, 2;, one can express the solution, regarded as a function of t, as 
the sum of two terms, one of which varies very slowly. The other term 
varies much more quickly and appears to be a combination of exponentials. 
The complex frequencies that appear in this latter combination of expo- 
nentials agree very closely with those obtained from a linearized version of 
the problem as described in Chapter 2 of part II of the above mentioned 
work. The error analysis in this report was based on the assumption that 
the functions to be integrated could be expressed in terms of these frequencies 
and the results seem to justify this assumption. (Cf. loc. cit., p. 121-123.) 

These frequencies are listed for various time values. The values are 
taken from tables on pages 72 and 73 of the report. dX. and sg are identically 
zero. The complex frequencies are listed in Table 2, the real frequencies 
in Table 3. The Table of § 10 may be used to compute the maximum A. In 


TABLE 1 
nm, M2 = —.35 + 5.6677 Mm, Ue = —.014 + .226682 
¥3,% = —.234 + 1.0647 Uz, Ue = —.00936 + .042567 
y= —2.9 ¥%= 0 “= —.116 um = 0 
TABLE 2 
t My, re As, M 
2.30 — .34965 + 5.664902 —.22970 + 1.065407 
2.70 — .34664 + 5.627102 — .22886 + 1.057081 
12.500 —.28158 + 5.071397 —.21667 + .950002 
30.020 —.19358 + 4.310712 —.18938 + .80357% 
50.020 —.11056 + 3.408467 —.14739 + .62710% 
t Ano, An Ar, Ais 
2.30 —.35142 + 5.668702 — .23858 + 1.062537 
2.70 —.35160 + 5.625017 —.23895 + 1.053272 
12.500 —.27892 + 5.062397 —.22938 + .946527 
30.020 —.19116 + 4.304927 —.20374 + .797981 
50.020 —.10608 + 3.406932 —.16836 + .62433% 
TABLE 3 
t As Au Ar As 
2.30 — 2.9009 — 2.8940 — .062634 +.0027338 
2.70 — 2.8980 — 2.8949 — .055760 +.0027604 
12.50 — 2.9165 — 2.9124 — .042666 +.0029659 
30.020 — 2.9496 — 2.9485 — .040503 +.0035029 
50.020 —2.9761 — 2.9761 — .044609 +.0042092 
TABLE 4 
t DY a 
0 2.97013 88889 2.90460 07614 
& — 5.50208 33333 — 5.19498 93728 
2 6.93194 44444 6.3 5668 29407 
3 — 5.06805 55556 —4.5 2952 65434 
4 1.99791 66667 1.74595 45079 
5 — .32986 11111 — .28272 22943 
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this case the smallest 4) obtainable from the table is associated with 


This yields 4) = .178. (The value of kp corresponding to Ay, Az or Ay is 
essentially the same. The values for A3, Ag = —.23 + 1.0654 is .85 and for 
the real A; = —2.9 is .239.) The experience with the problem seems to 


indicate that even h = .12 is too large when polynomials are used but it is 
conceivable that, with a proper choice of the »; even h = .08 might be 
effectively used with the exponential method. For comparison with the 
polynomial results let h = .04. 

The process of finding the a;’s is illustrated by considering the use of 6 
coefficients for the ‘Lark problem’’ and a time interval of 4 = .04. The 
following values of v; and u; were chosen. 

This will be referred to as the » = 6 case. In Table 4 the D; and the a 
appear so that the amount of change of coefficients needed is indicated. 


5 
Table 5 contains the yj, g(y;), > Ary7, hj, Tins(y; — ys). The values of 
r=0 


TABLE 5 3 
j i (93) 
1 — .01184 4296 — .22791 22417 .98457 97220 + .11180 859367 
2 — .01184 4296 + .22791 22417 .98457 97220 — .11180 859367 
3 .00848 9886 — .04294 72631 .99503 48169 + .02114 44856: 
4 .00848 9886 + .94294 72631 -99503 48169 — .02114 448563 
5 .12299 5872 .94417 91121 
6 0 J 
5 
> Ary 
j r=0 hy 
1 .98462 35536 + .11180 451322 — .00004 38316 + .00000 408043 
2 .98462 35536 — .11180 451327 — .00004 38316 — .00000 408043 
3 .99503 48145 + .02114 449167 .00000 00024 — .00000 000607 
4 .99503 48145 — .02114 449167 .00000 00024 + .00000 00060% 
5 .94417 81327 .00000 09794 
6 1. 0 
; 6 
TL (yy; — x) 
é=1 
3 
Jj 
1 — .00096 98667 — .00099 901727 
2 — .00096 98667 + .00099 901723 
3 .00002 01775 + .00001 153767 
4 .00002 01775 — .00001 153767 
5 .00012 89991 
6 — .00001 22777 


the p,,; and E;, ; are given in Tables 6 and 7. The values of K are given in 
Table 8. The a;’s are good to at least 7 figures after the decimal point. 
The case in which one is permitted only four coefficients may also be 
considered. In this case, one may choose »; and v2 to be —.292 + 3.3655 
corresponding to an average of the previously used ». The error formula 
shows that for d; and A:, the truncation error is reduced by 28 percent over 
the polynomial case. This relatively slight improvement is due to the need 
of approximating 6 complex numbers by 4 in a certain way. If one were 
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TABLE 6 
Pej 


N 


~ 


0 
2 
3 
os 
5 
0 
2 
3 
4 
5 


mA 


i 
j 
1 
j 
0 
1 
2 
3 
4 
5 


~ 


m. nm. mw. 
Aur whe AMP wWN-e Aur oN 
~ ~ ~ 


1 


1. 
-12813 1348 + .22791 22417 
.00234 70745 + .03190 21627: 
.00018 82909 +- .00091 27857: 
—_— 27920 + .00005 372514 


3 


:. 
10779 7166 + .04294 72634 
.05001 37854 + .00426 497591 
00682 36045 + .00211 174607 
Pats 43875 + .00027 51267: 


5 
1, 
— .00670 8820 
.05359 85956 
.00083 89789 
.00009 98219 
0 
TABLE 7 
Ej 
0 


1.13066 39214 + .26078 091452 
1.13066 39214 — .26078 091457 
1.16468 89434 + .04959 911162 
1.16468 89434 — .04959 91116: 
1.04782 85764 
1.17670 71658 


2 


10.77601 76025 + 1.46409 27198: 
10.77601 76025 — 1.46409 27198: 
10.80364 79567 + .27259 02517: 
10.80364 79567 — .27259 025171 
10.12138 18457 
10.86347 40121 


4 


5.12813 1348 + .22791 22417 
5.12813 1348 — .22791 22417 
5.10779 7166 + .04294 72634 
5.10779 7166 — .04294 72634 
4.99329 1180 
5.11628 7052 
TABLE 8 
K; 
.01982 50786 — 


2 


a. 
-12813 1348 — .22791 2241: 
00234 70745 — .03190 216274 
00018 82909 — .00091 27857: 
a 27920 — .00005 37251: 


4 


i. 
-10779 7166 — .04294 7263i 
05001 37854 — .00426 497592 
00682 36045 — .00211 174607 
oo 43875 — .00027 512673 


-11628 7052 

05277 34384 
.00743 13849 
-00020 30128 
00001 22777 


1 


5.51994 04053 + 1.00923 47486: 
5.51994 04053 — 1.00923 47486: 
5.59493 16167 + .189082 59843 
5.59493 16167 — .18909 259843 
5.13573 82865 
5.63853 43058 


3 


10.51487 24665 + .94355 11267: 
10.51487 24665 — .94355 112677 
10.48120 24494 + .17605 40279: 
10.48120 24494 — .17605 402793 
10.02676 33156 
10.51792 16464 


— ee 


.02462 812107 





01982 50786 + .02462 812107 
— .00005 18207 — .00027 53833% 
— .00005 18207 + .00027 53833% 
or 23010 


AMP aN &, 








aa sear A 


Padre &, 


PWN, 


TABLE 9 


D® 
2.29166 66667 

— 2.45833 33333 
1.54166 66667 

— 0.37500 00000 


@OnrK oO = 


TABLE 10 


Vi 
.00259 45860 — .13579 056607 
00259 45860 + .13579 056607 
~~ 58721 


3 
y Anyi 
r=0 


-.99107 63783 + .06666 376651 
.99107 63783 — .06666 37665% 
-94410 76360 
-00000 00000 


k=4 
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as 
2.23915 28992 
— 2.30866 86852 
1.39918 51376 
— 0.32966 93515 


g() 
.99119 12907 + .06668 754073 
-99119 12907 — .06668 754077 
.94417 91095 
.00000 00000 


hy 

00011 49124 + .00002 377424 
00011 49124 — .00002 37742¢ 
00007 14735 

00000 00000 


Tl (9; — y) 
k= 
kei 


00434 44938 + .00509 254534 
00434 44938 — .00509 254532 


an So wn oO 
— 
~~. n~. 


nm ~, 
mw ed 
~ ~~ 


-00405 09364 
— .00226 87585 


BON RR, 


TABLE 11 
Pri 
1.00000 00000 


0.12559 04581 + .13579 056607 
0.00031 91234 + .01670 167911 
0.00000 00000 


3 
.00000 00000 
0.00518 91720 
0.01844 58097 
0.00000 00000 
TABLE 12 
E,, 
0 +. 


1.12590 95815 + .15249 224514 
1.12590 95815 — .15249 224514 
1.02363 49817 
1.14953 78590 


2 


3.12559 04581 + .13579 056607 
3.12559 04581 — .13579 056604 
3.00518 91720 
3.12818 50441 


2 


1.00000 00000 
-12559 04581 — .13579 056607 
.00031 91234 — .01670 167914 
0.00000 00000 


1.00000 00000 
-12818 50441 
01908 40564 
00226 87585 


1 


3.25150 00396 + .28828 281114 
3.25150 00396 — .28828 281114 
3.02882 41537 
3.27545 41446 


— he pe ee 
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TABLE 13 
K; 
.01384 34770 — .01075 483487 


01384 34770 + .01075 48348: 
01764 36945 
-00000 00000 


PWN RR, 


permitted 6 coefficients, initially the truncation error practically vanishes. 
Even after a change of 10 seconds in the independent variable, the trunca- 
tion error has been reduced by 76 percent relative to the polynomial case. 
Tables 9 to 13 for m = 4 are analogous to Tables 4 to 8, respectively, 
for the n = 6 case. 
P. Brock 


Reeves Instrument Corp. F. J. Murray 
Columbia University 


The preparation of this paper was assisted by the Office of Naval Research. 


RECENT MATHEMATICAL TABLES 


992[C].—G. W. SPENCELEY, R. M. SPENCELEY & E. R. Epperson, Smith- 
sonian Logarithmic Tables to Base e and Base 10. Smithsonian Misc. 
Collections, v. 118, Washington, 1952, xii + 402 p., 14.6 X 22.9 cm. 
Price $4.50. 


This volume gives both common and natural logarithms of numbers of 
the forms 


n, 1+ 2-10, 1+ 2-10°, 


where m = 1(1)10*. 

The values are given to 23D. In the case of common logarithms the 
characteristics are omitted. P. xii contains 23D values of the natural log- 
arithms of 10*, = 1(1)10 as well as log e. 

The table is intended to be used with a calculating machine to find 
logarithms and antilogarithms to accuracy not exceeding 23D by the well- 
known factorization method. It was computed on desk calculators with the 
assistance of 4 students. The natural logarithms of. integers were built up 
from the table of WoLFRaAm. All common logarithms were found by multi- 
plication by log e. The work was carried to 28D and then rounded to 23D. 

These tables prove that it is still possible to produce a hand-set volume 
of over a million digits from a very small computing organization. The 
existence of large scale computing units that could calculate the present 
table in three days does not seem to daunt the authors. 

The FMR Index lists only one table comparable with the present one: 
the four-figure radix table of STEINHAUSER? to 21D which “contains many 
errors.” The need for the present table is certainly not as great as it was 
in 1880 with the advent of the modern automatic desk calculator. Having 
both natural and common logarithms is something of ‘a iuxury in the face 
of present day printing costs. The hand computer who has a need for many 
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very accurate logarithms will be grateful to the Smithsonian Institution for 
bringing out this new addition to their series of tables. 
D. H. L. 
1 The reprint used was that of Veca [MTAC, v. 4, p. 194]. 


2A. STEINHAUSER, Hilfstafeln zur pricisen Berechnung zwanzigstelliger Logarithmen- 
Vienna, 1880. 


993[D, E, K, L]|.—BAASMTC, Mathematical Tables, v. 1: Circular & Hyper- 
bolic Functions, Exponential & Sine & Cosine Integrals, Factorial Func- 
tion & Allied Functions, Hermitian Probability Functions. Third ed. 
Cambridge University Press, for the Royal Society, 1951, xl, 72 p. 
20.9 X 28 cm. 18 shillings. 


The first edition of these tables was published in 1931; in the review of 
the second edition, 1946 [MTAC, v. 2, p. 122-123] the list of the XVI 
tables (p. 1-72) was given. Most of the valuable Introduction of the first 
edition was omitted from the second edition, but it is now with suitable 
modification restored (p. viii-xxxvii) and elaborated (p. xxxviii-xl) for the 
third edition. Readers may be reminded that different sections of this Intro- 
duction had been prepared by A. J. THompson, L. J. Comrie, J. HENDERSON, 
A. Lopce & J. WisHart, J. O. Irwin, & R. A. FisHER. Thus the four-page 
Introduction of the second edition has now been dropped. 

On p. xxxviii the extension of Table II, Factorial Function, to 18D, and 
Constants, is reprinted from the second edition. The four new Constants 
now added are 16D values of 1/2, 2’, log x, In (2x)*. The remaining entirely 
new material of the Introduction includes: the values of the first 12 Ber- 
noullian numbers; Corrigenda to the first edition, 1931, the only serious 
error listed being in sin 47.6 [MTAC, v. 2, p. 135]. The serious error in 
cos 48.6 listed in 1952 for the second edition by S. H. Conn [MTAC, v. 6, 
p. 100] occurs also in the first and third editions. Then follow : Corrigenda 
in the first, second, and third editions; 9 unit errors in the final digits dis- 
covered by J. W. WRENCH, JR.; Corrigenda to the second edition; Addenda. 

This third edition, so much more satisfactory than the second, is indeed 
wholly admirable, except for that serious error (to which we referred), 
discovered too late for correction. 


R. C. ARCHIBALD 
Brown University 
Providence, R. I. 


994[D, R].—B. Goussinsky, Tables for Checking Traverse Computations with 
the aid of a Calculating Machine. [Palestine] 1952, Ministry of Labour, 
Survey of Israel, ii + 18 p., 14.6 X 23.5 cm. Price 400 pruta. 


The table gives 5D values of cos x — sin x for x = 0(1’)45°. Tables of 
proportional parts are intended for interpolating to 10”. 


995[F, G].—F. C. AuLuck, “(On some new types of partitions associated 
with generalized Ferrers graphs,” Camb. Phil. Soc., Proc., v. 47, 1951, 
p. 679-686. 
The author tabulates three special partition functions which have com- 
plicated verbal definitions. The first two may be defined in terms of the 
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number (m) of unrestricted partitions of n as follows: 


P(n) = p(n — 1) — p(n — 3) + p(n — 6) — p(n — 10) + --- 
Q(nm) = P(i)p(m — 1) + P(2)p(m — 2) + --- + P(m)p(0). 


Here 1, 3, 6, 10, --- are the triangular numbers. If we denote by p;(n) the 
number of partitions of m into precisely k positive parts, the third partition 
function may be defined by 


R(n) = Lbi(m)p(n — 2k? — k — 2m — 1), 


the sum extending over all positive integers (m,k). All three functions 
P, Q, R are tabulated for m = 1(1)20. 
D. H. L. 


996[F, G].—J. RrorDan, ‘‘The arithmetic of ménage numbers,” Duke Math. 

Jn., v. 19, 1952, p. 27-30. 

The general ménage number U,,, is the number of arrangements of n 
married couples seated at a circular table, gentlemen alternating with ladies, 
in which precisely r husbands are seated next to their own wives. The paper 
gives a small table of U,,., for r = 0(1)5 and k = 0(1)5. The heading of 
the table being a little obscure, it may be worth noting that, for example, 
Us,2 = 40. 

D. H. L. 


997(I].—H. E. Sauzer, “Formulas for numerical integration of first and 
second order differential equations in the complex plane,” Jn. Math. 
Phys., v. 29, 1950, p. 207-216. 


Let 
(1) Fa) = [" fis)ds; Ge) = f ” F(s)ds, 


where 2 = x + iy. If f(z) can be approximated by a polynomial in z of 
degree m — 1, then both F(z) and G(z) can be expressed as follows: 


F(a) & Fin = x ari f(t) + A’ 
(2) 
G(s;) = Gin = x biflm) +A +B. 


In the above A’, A, and B are constants of integration and the points z are 
fixed, convenient points in the complex plane at which f(z.) is known. The 
author tabulates the coefficients a,’ and 5’, for 7 = 1, 2, ---, , such that 
(2) coincides with (1) when f(z) is a polynomial of degree no higher than 
n — 1, for » = 2,3, ---,8. The points 2 were chosen over an L-shaped 
grid in the complex plane from the following considerations: 


(a) The configuration shall be such that it will be convenient to initiate 
integration in either the x- or y-direction. 
(b) The points 2, shall be as close together as possible. Thus for the 
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five-point formula, the choice of the grid was 
Zo, 20 + h, 2 + 2h, 2 + th, % + (1 + ah. 


This grid was selected in preference to one that is symmetric about the 
45°-ray even though the choice involved more tedious computations in the 
evaluation of the coefficients. 

(c) The configuration shall permit of simple translations in either the 
x- or the y-direction. 


Since, in stepwise integration of a differential equation, it is often neces- 
sary to start with some approximation to f(z) at a point where f(z) is not 
yet known, the author also provides the coefficients for the corresponding 
n-point extrapolation formulas for f(z), which are exact when f(z) is a poly- 
nomial in z of degree no higher than m — 1. 

With modern high-speed computing equipment, the evaluation of the 
constants of integration at each step of the calculation is manageable, and 
the formulas can be expected to give a better approximation, for a com- 
parable amount of work, than that which would result from breaking up 
f(z) into real and imaginary components, and then considering the integra- 
tion problem in real variables x and y. 

The paper occupies only 10 pages, but any one who has tried to calculate 
at least one set of the coefficients will appreciate the prodigious amount of 
work that this table represents. Because of the author’s reputation for 
accuracy in his many tables of coefficients, the present work can no doubt 
be used with complete confidence. Those interested in computing methods 
owe a debt of gratitude to the author for providing a well planned and useful 
new tool for numerical integration in the complex plane. 


GERTRUDE BLANCH 
National Bureau of Standards 


Los Angeles 24, Calif. 


998[K].—Nits BLomgvist, “On a measure of dependence between two 
random variables,” Annals Math. Stat., v. 21, 1950, p. 593-600. 


The author investigates the statistic, g’ = (m — m2)/(m +m), as a 
measure of dependence between two random variables. With the x — y 
plane divided into four quadrants by the axes x = ~ and y = §, m is the 
number of sample points in the ist and 3rd quadrants, while m, is the 
number in the 2nd and 4th quadrants. This statistic was suggested earlier 
by MosTELLER! as a “rough” but more readily computed measure of 
dependence than the Pearson correlation coefficient. He did not, however, 
fully investigate its properties. 

To facilitate tests of significance based on gq’ the author derives the 
necessary formulas and presents, in tabular form, values of P{|m, — k| > v} 
for »y = 0(1)15 and 2k for k = 2(1)25, where 2k is the largest even number 
contained in m (m = m, — m2). Separate tables are given for odd and even 
values of vy and k. Most entries are given to 3D. Those for smaller proba- 
bilities are carried to 4D. 

A. C. CoHEN, JR. 
University of Georgia 
Athens, Georgia 


1F, MosTE._er, “On some useful ‘inefficient’ statistics,” Annals Math. Stat., v. 17, 
1946, p. 377-408. 
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999[K].—Z. W. BrrnBaum & F. H. TinGey, “‘One-sided confidence contours 
for probability distribution functions,’ Annals Math. Stat., v. 22, 1951, 
p. 592-596. 


Let X be a random variable with the continuous probability distribution 
function, F(x) = Probability {X < x}. An ordered sample X; < X2 < 
< X, of X determines the empirical distribution function, F,(x) = 0, for 


x < Xa; Fale) = * for Xe Se < Xe b = 1,2, ---,2 — 1; F,(x) = 1 for 


Xn < x. The function F,(x) = min [F,(x) + ¢, 1], also determined by the 
sample, is called an upper confidence contour. Let P,(e) = Prob. { F(x) 
< Ft ,(x) for all x}. The authors derive an explicit expression for P,(e) and 
tabulate values of €,,. for a = .10, .05, .01, .001, and = 5, 8, 10, 20, 40, 50 
to at least 4D in Table 1. In Table 2 the authors present the values to 4D 
for the same values of ” and a as in Table 1 of &,4 = (2 In a)! obtained 
from the asymptotic results of N. Smrrnov.' For m = 50 the difference 
between €,,. and é,,2 is small. 


L. A. AROIAN 
Hughes Aircraft Company 


Culver City, California 


1N. Smirnov, “Sur les écarts ‘. la courbe de distribution empirique,” Rec. Math. [ Mat. 
Sbornik], n.s., v. 6, 1939, p. 3-26 


1000[(K].—E. G. CHAMBERS, Statistical Calculation for Beginners. Cambridge 
University Press, London and New York; second ed., 1952 (first ed., 
1940). x + 168 p. 14.9 K 22.2 cm. 12s. 6d. 


The title of this volume may be a bit deceptive, especially to readers of 
MTAC, since it is not concerned with computational methods as such, but 
is rather another introductory and non-mathematical book on statistical 
methods for “students of the biological sciences, especially psychology.” 
The details of the calculation of elementary statistical quantities and stand- 
ard significance tests (there are some for use with qualitative data which 
not all statisticians would regard as standard) are explained and illustrated 
and there are further exercises with answers for the reader to do. The 
student is not assumed to have access to a desk calculator. If this 
review were addressed to a statistical audience certain inaccuracies of 
statement should be mentioned though these are not numerous or really 
serious. Besides standard statistical tables, usually in a shortened form, 
there are some auxiliary tables for assistance in computation, namely of 
[nyn2(m + m:— 2)/(m + nz) |§ for ny and nz = 10(1)50 to 2D; N(N? bai 1)/6 
for N = 10(1)69; n(m — 1)(2n + 5) for m = 2(1)60; ¢#(t — 1)(¢ — 2) for 
t = 3(1)50; mn(m — 1)(m — 1)/4 for m = 2(1)6 and m = 3(1)15; 


(5mm — 1)/(m — 2)? for m = 3(1)6 and m = 3(1)15; (5)(7) (m — 3)/ 


2(m — 2) for m = 4(1)6 and m = 2(1)15; and m?(n* — n)/12 for m = 3(1)6 
and m = 3(1)15. 
a Tae Se 
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1001[K].—ArTHUR LINDER, Statistische Methoden fiir Naturwissenschafter, 
Mediziner und Ingenieure, Verlag Backhauser, Basel. Second ed., 1951 
(first ed., 1945). 238 p. 17.5 K 24.4 cm. 30 Swiss francs. 

This book filled a need for a modern introductory book in mathematical 
statistics, including derivations, in German. Besides standard statistical 
tables it includes a table for the direct evaluation of the significance levels 
of regression coefficients. 5, 1, and .1 percentage points for 1(1)6 independent 
variables and for 1(1)30, 40, 60, 120 degrees of freedom are given to 4S. 


Cac. 


1002 [K].—Sice1t1 Moricvut1, ‘Extremal properties of extreme value distri- 

butions,”’ Annals Math. Stat., v. 22, 1951, p. 523-536. 

Schwarz’s inequality is used to obtain an upper bound for the expecta- 
tion, a lower bound for the coefficient of variation, and lower and upper 
bounds for the standard deviation of the extreme value taken from a sym- 
metrical distribution with zero expectation and unit standard deviation. 
The bounds, given as functions of the sample size m, are actually reached 
for specific distributions depending upon m. The upper bound for the expec- 
tation has been already derived by PLACKETT,' by a more complicated 
method based on the calculus of. variations. 

Table 1 gives 4D and 5D values for the upper bound of the expectation 
of the extreme value, the (known) expectation of extreme values for the 
normal distribution and those for a rectangular distribution for m = 2(1)20. 
Table 2 gives to 4D and 5D the lower bound of the coefficient of variation 
of the extreme value, the (known) values for the normal and those for the 
rectangular distribution for nm = 2(1)6. Seven graphs show the upper bounds 
for the expectation of the extreme value for = 20, the lower bounds for 
the coefficient of variation for » = 10, the upper and lower bounds for the 
standard deviation of extreme values for » = 10 and the probability and 
density functions where the bounds are reached. The proofs are unusually 
simple and clear. 

E. J. GUMBEL 
New School for Social Research 
New York, N. Y. 


1R. L. PLackett, “Limits of the ratio of mean range to standard deviation,” Bio- 
metrika, v. 34, 1947, p. 120-122. 


1003[K, L].—P. Natu, “Confluent hypergeometric function,” Sankhyd, v. 
11, 1951, p. 153-166. 


The present paper contains 7S tables of 
M(a, y,x) = DT) (a + v)x*/[r'P(@)ly + 1)J 
r=0 


for y = 3, a = 1(1)40 and y = 4, a = 1(1)50; in both cases x = .02(.02)- 
.1(.1)1(1)10(10)50, 100, 200. The computation of the tables is described in 
the introduction and the author states that as he worked throughout with 
10S and as the 40 times repeated use of the recurrence relation produced 
errors in at most the last two figures, the seven significant figures given in 
the table “‘can be taken to be all correct.” 
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It is pointed out that in statistics special forms of M(a, y, x) occur in 
the distribution functions of the F statistic and D? statistic and in the error 
function and the incomplete function. 


I. R. SAVAGE 
NBSSEL 


1004[K].—PauL R. Riper, ‘The distribution of the quotient of ranges in 

samples from a rectangular population,’’ Amer. Stat. Assoc., Jn., v. 46, 

1951, p. 502-507. 

Two independent random samples of size m and m from a continuous 
rectangular population have ranges R; and R:, respectively. The density 
function of the quotient, « = R,/R: is given and upper percentage points 
for sample sizes not greater than 10 are computed. Tables 2, 3, and 4 give, 
respectively, the upper 10, 5, and 1 percent points of the distribution of u 
to 3S for m and n = 2(1)10. 


LEo Katz 
Michigan State College 


East Lansing, Mich. 


1005[K].—W. L. STEvEns, “‘Asymptotic regression,’’ Biometrics, v. 7, 1951, 

p. 247-267. 

In fitting a regression curve of the form y = a + p* by least squares, 
estimates a, b and r are obtained of the three parameters. A covariance 
matrix can be written in terms of 5 and of six functions Fas, Fas, Far, For, 
Fy, F,, of r alone, namely, of the components of the inverse of 


n zr? 2xr?-! 
=r? zr 2xr*=-! 
Dxr2-! Ter) Yx2z22-2 7, 


For the case where the x values are spaced equally over m values (coded 
as 0,1, ---,#— 1) values of the functions are tabled to 5D: m = 5, 6, 
r = .25(.01).65(.005).70, and m = 7, r = .30(.01).70(.005).75. 

The tables facilitate the numerical problem of approximating values of 
a, b, r which satisfy the least squares principle. Examples are given. 

F. J. MASSEY 

University of Oregon 
Eugene, Ore. 


1006[K].—MaryjoriE Tuomas, ‘‘Some tests for randomness in plant popu- 

lations,” Biometrika, v. 38, 1951, p. 102-111. 

Power functions are obtained for three tests designed to detect departure 
from randomness in the distribution of a plant species over an area. The 
area is divided into N quadrats, in each of which the probability, P;, of 
observing & plants is 
k e-™m’ e(ry)** 
et ar Be | 
on the assumption that plants are distributed in clusters, the number of 
clusters per quadrat having a Poisson distribution with mean m and the 
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number of plants in excess of one per cluster having a Poisson distribution 
with mean \. The hypothesis of randomness (A = 0) is to be tested against 
alternatives, \ > 0. 

The three tests compared are (i) STEVENS” test based on m, the number 
of quadrats containing no plants, (ii) s = }>(x; — Z)?/Z, where x; is the 


“ a eee m 
number of plants in the ith quadrat and (iii) ] = — In no In | . 


the maximum likelihood estimate of \ based on observed frequencies of 
quadrats with none and one plant. As might be expected, z is the preferred 
statistic. 

Table 1 gives the power of Stevens’ test for a = .03 to 3D for A = 0(.2)1.4 
when N = 25 exactly for two artificial situations in which the total number 
of plants observed, S, is taken to be 5 and 10, respectively, and approxi- 
mately for a more realistic situation where M = m(1 + A) = .4. It should 
be noted that the approximate powers were obtained by experimental 
sampling, in each case involving 200 random samples of 25. 

Table 2 gives the power of z for a = .05 to 3D for A = 0(.1).6; N = 100, 
200; M = 1, 2, 3. Table 3 gives the power of / for a = .05 to 3D for 
A = 0(.1)1.4; N = 100, 200, 300; M = 1, 2, 3. 

Leo Katz 
Michigan State College 
East Lansing, Mich. 


1W. L. Srevens, “Significance of grouping,’’ Annals of Eugenics, v. 8, 1937, p. 57-69. 


1007[K].—H. R. VAN DER VAaRT, “‘Gebruiksaanwijzing voor de toets van 
Wilcoxon,’’ Mathematisch Centrum, Statistische Afdeling, Rapport S 32 
(M 4), 1950. 

Let X and Y be independent random variables with continuous cumu- 
lative distribution functions F(x) and G(y). Let m, ---,x, be a set of a 
independent observations on X and 4, ---, ¥ma set of m independent obser- 
vations on Y. By u we denote the number of pairs (7, j) with x; > y;. The 
quantity u is used in WILCoxon’s test! for the hypothesis H that X and Y 
are the same random variables (F(x) = G(y)) against the alternative H that 
F(x) > G(x) for all x (Y is stochastically larger than X). Whenever u = tu 
the hypothesis H is rejected, u being chosen such that under H the proba- 
bility P(u S um) is about equal to the desired level of significance a. This 
probability depends only on m%, m, and m; it will be denoted here by 
H(uo, n,m). It satisfies H(u, m,n) = H(u, n,m), H(mn — u, m,n) = 1 
— H(u, m,n) and H(u, m,n) = 0 for u < 0. Hence in tabulating this func- 
tion only integral values of m, n, and u need to be considered with 1 < m <2 
andOQ<u< = . Complete tables of P(u = u) = H(u, m,n) have been 
given by MANN & WHITNEY? for m, m = 1(1)8. The present report gives 
complete tables for m, m = 1(1)10 to 7D. 

Wilcoxon's test is important in comparing different treatments (with 
“effects” X and Y) when, except for continuity, nothing is known about 
F(x) and G(y). The research worker without any statistical background can 
find in this report a very detailed explanation of this test and the relatively 
easy computations involved. The author also points out how to obtain 
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P(u S u) whenever m and n are large, using the normal distribution as an 
approximation. 
J. H. B. KEMPERMAN 
Purdue University 
Lafayette, Ind. 
1 FRANK WILCOxON, “Individual comparisons by ranking methods,” Biometrics Bull., 
v. 1, 1945, p. 80-83. 


2H. B. Mann & D. R. Wuitney, “On a test of whether one of two random variables 
is stochastically larger than the other,” Ann. Math. Stat., v. 18, 1947, p. 50-60. 


1008[K].—J. E. Watsu, ‘‘Some bounded significance level properties of the 
equal-tail sign test,’”” Annals Math. Stat., v. 22, 1951, p. 408-417. 


The tables exhibit coefficients for confidence intervals based upon a 
simple equa!-tail procedure (essentially a sign test) for testing hypotheses 
about the population median. Significance levels are determined by confi- 
dence coefficients which depend upon the population conditions. 

Two tests are used. Test 1, the simpler, assumes two conditions on the 
populations: (i) the populations have a common median value, (ii) no 
population has a discrete amount of probability concentrated at the median. 
Table 1 shows significance levels for confidence intervals for samples of a 
assuming conditions (i) and (ii) for m = 4(1)15 to 4D. 

When conditions (i) and (ii) are not satisfied, the null hypothesis is 
stated in terms of some function of a set of median values contained in the 
set of 50% points common to all of the populations. Test 2 applies here and 
this reduces to Test 1 when the set becomes a single value, thus satisfying 
condition (i). Test 2 is not exact but rather gives upper and lower bounds 
to the significance level. These are aiso given in Table 1. If the populations 
are assumed to be continuous at the median, the lower bounds of Test 2 
are considerably improved. These lower bounds are given in Table 2 for 
nm = 4(1)15 to 4D. 

T. A. BICKERSTAFF 
University of Mississippi 
University, Miss. 


1009[L].—V. A. Ditkin & P. I. Kuznetsov, Spravochnik po operafsionnomu 
Ischislenitu Osnovy Teorii i Tablitsy Formul. [Handbook of operational 
calculus. Fundamentals of the theory and tables of formulas.] Moscow 
and Leningrad, Gostekhizdat, 1951, 256 p. 14.6 X 22.2 cm. Price 7.30 
roubles, bound. On pages 98-105 are definitions of higher functions. 

In the Heaviside operational calculus we can use the Laplace transform. 

The transform of a function f(t) is 


Fe) = p f° sera, 


provided f(t) is of exponential order at infinity, and then has the inversion 
formula 


at+io 
HO) = 35, [eR eran. 
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On pages 106-111 are general formulae listing F(p) and corresponding f(é). 
Then on pages 112-255 are tables of transforms F(p), f(t) arranged under 
the following 15 function headings : Rational; Irrational; Exponential; Trigo- 
nometric, and hyperbolic; Logarithmic, inverse hyperbolic, and inverse 
trigonometric; Gamma and related; Integral; Confluent hypergeometric; 
Bessel; Legendre; Elliptic; Theta; Mathieu; Hypergeometric series; various. 
In all there are 1354 transform pairs. 

There is no reference to the elaborate similar lists in (a) J. Cossar & 
A. ErvE yt, Dictionary of Laplace Transforms (1944-1946; see MTAC, v. 1, 
p. 424-425; v. 2, p. 76, 215-216), or in (b) N. W. McLacuian, P. Hum- 
BERT & L. PoLt, Supplément au formulaire pour le calcul symbolique (1950) 
although the Second edition of the original Formulaire is mentioned. 

This work was published in an edition of 10,000 copies. 


R. C. ARCHIBALD 
Brown University 


Providence, R. I. 


1010(L].—Otto EMERSLEBEN, ‘‘Numerische Werte des Fehlerintegrals fiir 
nz,” Zeit. angew. Math. Mech., v. 31, 1951; p. 393-394. 


F(x) = = f e~* dt. On p. 393 is a table of F( nx), for mn = [1(1)11; 


15D], for m = 1, 17D. 


R. C. ARCHIBALD 
Brown University 
Providence, R. I. 


1011(L].—G. B. Hacen, “Uber iterierte Integration von Bessel-Funk- 
tionen,” Zeit. angew. Math. Mech., v. 32, 1952, p. 27-30. 


On p. 29 is a table of f * Jo(t) dt, x = [0(.2)15(1)24, ©; 4D). 


0 
This integral has been tabulated to 8D by A. N. Lowan & M. ABRAmo- 
witz for x = 0(.01)10 [see MTAC, v. 1, p. 154]. From this table it may 
be verified that in Hagen’s table there are 13 fourth-place unit errors for 
x = .8, 1.8, 2.4, 3.6, 4.2, 4.4, 5.2, 5.6, 6.4, 6.6, 7.4, 7.6, 8.6, 8.8; fourth-place 
2-unit errors for x = 4.0, 4.8, 5.4; fourth-place 3-unit errors for x = 9.2, 9.8; 
fourth-place 60-unit error at x = 3.2. 


R. C. ARCHIBALD 
Brown University 
Providence, R. I. 


1012[L].—J. M. HammersLey, “On a certain type of integral associated 
with circular cylinders,’ Roy. Soc. London, Proc., v. 210A, 1951, p. 
98-110. 


Table 1, p. 104. 8D table of 


Bo) = 9 R(—8, 425) — 1 4 Ct ae} 


for c = 1(1)10. 
Table 2, p. 106, 107. 7D table of 


g(t) = 16¢{1 — EF(—4, 4; 2; #)} = 16% {1 — F(—4,452;E>)}, 
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h(é) = — 2 + 4r€{(1 — 4%) cos & + &(2% + 1)(1 — #)4} 


for = 0(.01)1, with second central differences (some modified), and in the 
case of g() also a few fourth differences. 

Table 3, p. 108, 109. 7D table of g(é) for = 1(.01)2(.05)3(.25)10, with 
second central differences (some modified) and a few fourth differences. 

k(c) has been computed from its expansion in descending powers of c, 
and g(é) from its expression in terms of complete elliptic integrals. The 
elliptic integrals were partly taken from existing tables,! and partly they 
were computed by the author. 

A. E. 


1 A. FLETCHER, “A table of elliptic integrals,” Phil. Mag., s.7., v. 30, 1940, p. 516-519. 


1013[L].—L. Howartn, ‘“The boundary layer in three-dimensional flow,” 

Phil. Mag., s. 7, v. 42, 1951, p. 239-243, 1433-1440. 

Part I of this work considers the problem of calculating the boundary 
layer for a fluid of small viscosity (large Reynolds number) flowing over a 
surface of quite general curvature. In Part II the work is extended to 
include the boundary layer near the stagnation point of a three-dimensional 
body where the stagnation point is a regular point of the body. 

Because it is a regular point of the body, a pair of coordinate directions 
x, y may be found such that, in the immediate neighborhood of the stagna- 
tion point, the flow at each axis is only in the directicn of that axis; that 
is, x, y are principal axes. The free stream velocity components U, V in 
the directions x, y may then be given 


U = ax; V = by (a, b constants) 


so that the flow is characterized by the ratio c = b/a. Let the corresponding 
velocity distributions in the boundary layer be 


u =axf'(z); v= byg’(s), 
where 


s = E(a/r)}, 


v being the coefficient of kinematic viscosity and £ the distance normal from 
the surface of the body. Howarth shows that 


(f')? — ff" — caf” = 1+" 
1 1 
(g")? beds ge” sk — =1+ -_ 
with boundary conditions 
f(0) = f'(0) = g(0) = g’(0) = 0 


f'=2'-1 as 3 &, 


Expanding in powers of c, 


f=frotcht+ eho, 
£ = Bot chi t+ cg2. 





a, «ss = SO. 4 OU 


Ls a a 7 


a= A ~*~ AO HS 





RECENT MATHEMATICAL TABLES 161 


Howarth solves successively for fo’(z), go’(z), fr’(z), gi’(z), fe’(z), ge’(2) by 
numerical methods. Table I gives these results to 3 decimals for 0.0 < z < 4.0 
at intervals of 0.1. The values of f’(z) and g’(z) are then calculated to 3 
decimals for ¢ = 0.25, 0.50, 0.75, 1.00 and given in Table II for 0.0 < z < 3.0 
for intervals of 0.1. The values of f’’(0) and g’’(0) are required for calcu- 
lating skin friction and are given in Table III to 3 decimals for c = 0, 0.25, 
0.50, 0.75, 1.00. Finally in Table IV are given values to 3 decimals of f(z) 
and g(z) for 0.0 < z < 4.0 at intervals of 0.5 in z. 


F. E. MARBLE 
California Institute of Technology 
Pasadena 


1014[L].—Davip MippLeton & VirGIniA JOHNSON, “A tabulation of 
selected confluent hypergeometric functions,’ Cruft Lab. Tech. Report 
no. 140, Jan. 5, 1952. 


The tables give 


\Fi(a;¢; —p) = Ere + n)P(0(—p)*/{T(@)T(c + n)n!}, 


generally to 5S, for p = 0(.25)2(.5)10(10)100 and c = 1, 2,a@ = —.5(1)8.5; 
c = 3, 4,5, a = —.5(1)9.5; c =.6, 7, a = —.5(1)10.5; c = 8, 9,a = —.5 
(1)11.5; ¢c = 10, a = —.5(1)12.5. A brief account of some properties of the 


confluent hypergeometric function and of some of its applications in noise 
problems is included as is a description of the computation. The bibliog- 
raphy of tables is hardly adequate. 

A. E. 


1015[L].—E. O. PowE t, “‘A table of the generalized Riemann zeta function 
in a particular case,’”’ Quart. Jn. Mech. Appl. Math., v. 5, 1952, p. 
116-123. 


Since the Hurwitz zeta function ¢(s, a) satisfies 
c(s, a) pe od t(s,a + 1) =a", 


a study of the difference equation f(a) — f(a + 1) = a led the author to 
a study of ¢(4, a). By means of the equation 
oo =(2n)! 
(1) (3,1+a) = Da s(n + 4)(-@)", = [a| <1, 
n=0 2 (n!) 
the calculation, at least for small a, is reduced to Riemann zeta functions 
which have been tabulated by Gram.! 
The multiplication theorem 


rd, Ma) = art E e(4,0+ 77) 


then permits calculation for larger values of a, and the difference equation 
serves as a check. A ten-decimal table of ¢(},a) is thus obtained for 
a = 1.00(.01)2.00(.02)5.00(.05)10.00 along with modified second central 
differences. 

The author also refers to an unpublished five-decimal table by R. HENs- 
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MAN? which gives {(s, a) for s = —10(.1)0, a = 0(.1)2, and (s — 1)¢(s, a) 
for s = 0(.1)1, a = 0(.1)2, obtained from a modification of equation (1). 


; T. M. Apostot 
California Institute of Technology 


Pasadena 4, California 


1J. P. Gram, ‘‘Tafeln fiir die Riemannsche Zetafunktion,’’ K. Danske Vidensk. Selskab, 
Skrifter, s. 8, v. 10, 1925, p. 311-325. 


?R. HENSMAN, Tables of the ae Riemann zeta function, Report No. T 211, Tele- 


— Research Establishment, Ministry of Supply, Great Malvern, Worc. 


1016[L].—L. WEINBERG, “Solutions of some partial differential equations 
(with tables),’’ Franklin Institute, Jn., v. 252, 1951, p. 43-62. 
The author considers the following differential equations: 
a) Laplace’s equation V*¢ = 0 
b) Poisson’s equation V7@ = — K 
c) Diffusion equation V*¢ = h-*d¢/dt 
d) Wave equation V*¢ = c*d¢/af? 
e) Helmholtz’s equation V*6 + 7 = 0, y* = constant. 


It is pointed out that a)-d) are special cases of Helmholtz’s equation. This 
is solved by the well known method of the separation variables. & is repre- 
sented as a product of three functions U', U?, U*, where each of those func- 
tions depends only on one variable. In the tables (p. 50-62), the expressions 
for V, V*&, div E, curl E and the three product functions U', U?, U*, are 
given (the latter also for the case of Laplace’s equation, y = 0) for the 
following system of coordinates: rectangular, circular cylinder, elliptic 
cylinder, parabolic cylinder, spherical, prolate spheroidal, oblate spheroidal, 
and parabolic. 


F. OBERHETTINGER 
American Univ. 


Washington, D. C. 


1017(L].—J. E. WiLxins, Jr., Some Miscellaneous Mathematical Problems. 
U. S. Atomic Energy Commission, Tech. Info. Service NYO-641, Oak 
Ridge, Tenn., Dec. 19, 1951, 21 mimeographed pages. 

The table of Laguerre polynomials described in UMT 135 [MTAC, v. 5, 

p. 232] is reproduced here on p. 15-21. 


1018[Q].—RoseErt J. Davis, Table of the Secant of the Zenith Distance for 
the Latitude of Oak Ridge, Massachusetts (42°5 N). Harvard Observatory 
Mimeograms, Series IV, No. 1, Cambridge, Mass., 1951. vi, 35 p. 
This table is for astronomical use, in correcting photoelectric measure- 

ments of stellar magnitude for atmospheric extinction. It gives to 3D 


sec z = (sin ¢ sin 6 + cos ¢ cos 6 cos t)-! 


for declination 6 = —28°(1°)90°, hour angle ¢ = 0°(4™)12" and latitude 
@ = +42°5. Auxiliary tables allow extension to any latitude +40° to +45°, 
and supply corrections arising from curvature of the earth’s atmosphere. 


JosEPH ASHBROOK 
Yale University 


New Haven, Conn. 
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MATHEMATICAL TABLES—ERRATA 
In this issue references have been made to errata in RMT 993 and 1011. 


210.—D. L. ARENBERG & D. Levin, Table of Fresnel Integrals and Derived 

Functions. Boston, 1948. 

In MTAC, v. 3, p. 479, Editorial Note, lines 8-10, I stated, ‘“‘The ArEn- 
BERG & LEvIN tables of C(u) and S(u), for u < 8.5, are identical with those 
given in JAHNKE & Empe, Tables of Functions, 1945, p. 34.”’ After “1945” 
I should have added, “‘p. 34, except in two cases.” For S(1.0) A. & L. had 
.4883 instead of .4383; for S(7.3) they had .3189 instead of .5189. 


R. C. ARCHIBALD 
Brown University 


Providence, R. I. 


211.—P. K. Bosr, “On recursion formulae, tables and Bessel function 
populations associated with the distribution of classical D*-statistic. 
Part I. On the construction of incomplete probability integral tables of 
the classical D*-statistic,’”’ Sankhyd, v. 8, 1948, p. 235-248. [MTAC, 
v. 3, p. 415.] 

The author lists in Sankhyd, v. 11, 1951, p. 96 a number of errors in the 
tables appearing in the above paper. 
.. oo 


212.—J. P. Kur, L. PoLetti & R. J. Porter, Liste des nombres premiers 
du onziéme million. |MTAC, v. 6, p. 81.] 
Page 23, Col. 5, Line 4 up. 
For 51 Read 511 
The number 10902511 is prime, so also is 10902517 which follows it. 


J. C. P. MILLER 
Univ. Math. Lab. 


Cambridge, England 
EpitoriaAL Note: This erratum occurs in one of our two copies but not in the other. 
D. H. L. 


213.—Rueticus-Pitiscus. A few months before he died Comrie sent to me 
a number of sheets headed “Children of Rheticus.’’ [Rheticus was born in 
1514 and died in 1576; see MTAC, v. 3, p. 552-562. ] From these sheets I 
abstracted certain tabular errata, checked every entry with the appropriate 
volume, and edited the whole in the form given below. The first three errata 
are in the 15D table of sines and cosines in Thesaurus Mathematicus (1613) 
of Rheticus, edited and published by BARTHOLOMAUs PitTIscus (1561-1613); 
see MTAC, v. 3, p. 390-398. The remaining errata are in the 10D table of 
tangents and cotangents in Opus Palatinum de Triangulis (1596) of Rheticus. 
For each angle in question the complete Rheticus value is given. Then in 
a new column are values derived from ANDOYER’s tables.' In the case of the 
first three values the corresponding Andoyer 15D values are transcribed. 
Corresponding to the 10D Rheticus tables, values got by rounding off 
Andoyer’s 15D tables to 10D values were found but only those Andoyer 
pentads corresponding to erroneous Rheticus pentads are indicated. Thus 









164 MATHEMATICAL TABLES—ERRATA 


it is seen that there are 164 erroneous entries, but 53 of these are simply 
unit errors in the tenth decimal place; see MTAC, v. 2, p. 285, lines 33-36. 

As explained in my Pitiscus article, the values given by Rheticus for 
tangents and cotangents less than 7° were especially faulty; hence in 1607 
Pitiscus published a new table for these. In this range we note below that 
there are 15 errors. On comparing the corresponding 1607 tabular values 
with those of Andoyer, I found that there was exact agreement in 12 cases, 
and that in the remaining three cases there was agreement to the extent 
that the tenth decimal place was one unit less in Pitiscus for 2° 57’ 50”’ and 
6° 45’ 20’’, but one unit more for 5° 26’ 40”. In the first of these three cases 
Pitiscus made no change in Rheticus; here was the first of the 53 cases of 
unit errors after rounding off Andoyer. 

So far as I am aware there exist only three original copies of the Pitiscus 
table of 1607; the Brown University film copy of Mr. W. D. MorGan’s 
original was the one I used in checking 15 values. For checking the 161 
values of the Opus Palatinum | used the excellent copy acquired in 1950 for 
the fine collection of old science books belonging to Mr. ALBERT E. LOwNEsS 
of Providence, R. I. The collation of this copy agrees in almost every par- 
ticular with that given in MTAC, v. 3, p. 556. 

R. C. ARCHIBALD 
Brown University 
Providence, R. I. 


1H. ANDOYER, Nouvelles Tables Trigonométriques Fondamentales, v. 1, 2, Paris, 1915, 
1916. 


Angle Function Rheticus Andoyer 
8 32 50 cos 0.98889 27054 25336 0.98889 37054 25337 
14 18 30 sin 0.24713 99747 34879 0.24713 99477 34878 
29 15 10 sin 0.48866 25407 22794 0.48866 35407 22794 
0 13 10 cot 261.09341 59146 47747 
0 44 10 cot 77.83149 34719 33743 
0 46 50 cot 73.39930 47413 44952 
1 20 00 cot 42.96407 72888 73358 
a tan 0.05177 58111 58112 
3 02 10 cot 18.85377 03160 02992 
5 26 40 cot 10.49202 10301 10272 
6 00 30 cot 9.50107 13501 13494 
6 24 30 cot 8.90351 05550 05496 
6 30 40 cot 8.76178 03398 03376 
6 38 10 cot 8.59529 72506 72473 
6 45 20 cot 8.44194 44514 44490 
6 47 00 cot 8.40705 11059 15030 
6 54 50 cot 8.24679 29153 28169 
6 59 30 cot 8.15415 08503 08485 
7 00 OO cot 8.14434 64279 64280 
7 02 10 tan 0.12342 42700 43701 
7 #03 20 cot 8.07957 25619 25482 
7 04 10 cot 8.06353 73496 73512 
7 06 SO cot 8.01264 33506 33469 
7 07 30 cot 8.00001 85504 85473 
7 14 S50 cot 7.86368 75505 75487 
cS, cot 7.74661 64501 64496 
7 39 30 cot 7.43688 70721 70701 
7 45 40 cot 7.33721 20493 20510 
7 52 10 cot 7.23495 06502 06488 
8 04 50 tan 0.14197 48487 48488 
8 05 10 cot 7.03859 70500 70508 
8 31 50 tan 0.14999 62513 62574 
8 


cot 6.59057 97491 97501 
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aly Angle Function Rheticus Andoyer 
; ° , ” 
om 8 49 30 cot 6.44101 64498 64507 
or 8 55 20 cot 6.36970 23503 23494 
907 8 55 30 cot 6.36768 74494 74505 
bat 9 15 30 cot 6.13460 62503 62499 
| 9 33 40 cot 5.93685 80614 80602 
ues 9 40 20 cot 5.86736 72518 72496 
es, 9 51 40 tan 0.17382 86080 86081 

t 10 01 50 cot 5.68364 92518 92500 
en 10 04 30 cot 5.62819 08495 08509 
ind 10 05 10 cot 5.62186 09485 09476 
ses 10 11 20 cot 5.56395 70500 70503 
of 10 24 30 cot 5.44411 18058 18053 
} O 10 30 50 cot 5.38422 74506 .38822 74493 

10 48 30 tan 0.19091 09931 09431 

— 10 49 30 tan 0.19121 24500 24501 
7 11 22 10 cot 4.97313 38625 38415 
N's 11 30 50 cot 4.90906 56814 56310 
161 1i 36 40 cot 4.86682 83675 83659 
f 11 39 00 tan 0.20618 01817 01187 
or 11 59 30 cot 4.70799 70513 70510 
TES 12 13 50 cot 4.61327 08501 08495 
age 12 23 50 cot 4.54931 25660 24661 
12 32 30 cot 4.49523 56502 56498 

12 39 40 tan 0.22464 66717 66317 

12 47 40 tan 0.22709 24603 24604 

12 55 30 cot _ 4.35749 13782 13770 

13 49 40 cot 4.06276 67383 67781 

13 53 10 cot 4.04501 70041 71040 

915 13 56 10 cot 4.02991 90504 90496 
’ 14 01 00 cot 4.00581 64036 65036 
14 08 20 cot 3.96976 11501 11498 

14 09 30 tan 0.25226 52596 52597 

14 31 30 cot 3.85976 46357 46758 

14 32 40 cot 3.85437 64095 65094 

$337 14 33 30 tan 0.25970 10737 40738 
1878 14 47 50 tan 0.26415 96316 95317 
1794 15 12 10 tan 0.27174 60475 60476 
15 13 30 cot 3.67417 43503 43497 

15 26 50 cot 3.61882 47481 47482 

15 30 50 cot 3.60249 21270 21965 

15 45 50 tan 0.28229 08741 08742 

16 00 50 tan 0.28700 77722 77422 

16 15 40 cot 3.42836 96577 96575 

16 23 50 tan 0.29426 33251 33252 

17 09 00 cot 3.24048 50395 60395 

17 23 00 cot 3.19425 97618 97619 

17 24 SO tan 0.31364 72372 72373 

18 28 10 cot 2.99399 02821 02323 

18 28 30 cot 2.99302 43787 43782 

18 32 00 cot 2.98291 66502 66505 

18 37 10 cot 2.96810 74167 74170 

18 46 40 cot 2.94121 94465 94965 

18 47 30 cot . 2.93888 17496 17500 

19 05 00 tan 0.34595 53113 53114 

19 09 40 cot 2.87789 68129 68130 

19 17 40 cot ~ 2.85643 95071 96070 

19 25 30 cot 2.83570 40500 40504 

19 27 00 cot 2.83176 39504 39499 

19 36 30 cot 2.80703 42717 42716 

19 46 10 tan 0.35941 98500 98501 

19 52 00 cot 2.76749 90375 90376 

20 04 40 cot 2.73591 59501 59498 

20 22 30 tan 0.37140 00666 00667 

20 40 10 cot 2.65069 74995 74998 

20 40 50 cot 2.64914 18153 18154 

20 41 50 cot 2.64681 12839 12834 

20 45 10 tan 0.37892 15500 15501 
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Function 


tan 
cot 
cot 
cot 
tan 
cot 
tan 
cot 
cot 
tan 
tan 
cot 
cot 
cot 
cot 
cot 
cot 
tan 
cot 
cot 
tan 
tan 
cot 
cot 
cot 
cot 
cot 
tan 
tan 
cot 
tan 
cot 
cot 
tan 
tan 
cot 
tan 
cot 
cot 
tan 
cot 
cot 
cot 
cot 
cot 
cot 
tan 
tan 
tan 
tan 
tan 
tan 
cot 
cot 
cot 
cot 
tan 
cot 
cot 
cot 
cot 
tan 
tan 
cot 


MATHEMATICAL TABLES—ERRATA 


Rheticus 


0.38403 09198 
2.59081 25803 
2.57445 12321 
2.55753 72314 
0.39738 68503 
2.47785 33464 
0.40509 81362 
2.37858 31847 
2.35394 83407 
0.42693 72163 
0.42951 90500 
2.29497 26583 
2.26343 88926 
2.21246 61499 
2.17138 23500 
2.16806 17501 
2.16172 06216 
0.47275 71499 
2.09863 79499 
2.06732 02054 
0.49664 60462 
0.50155 12499 
1.91961 85729 
1.88094 64517 
1.83709 54009 
1.82338 81499 
1.76051 70444 
0.57949 54234 
0.59149 37910 
1.68437 79623 
0.61353 46886 
1.62600 66589 
1.59053 70349 
0.62892 14878 
0.63995 66167 
1.55252 27347 
0.66684 70537 
1.45470 68333 
1.42800 16525 
0.71439 08583 
1.36230 37483 
1.35719 33743 
1.31070 91949 
1.29813 58518 
1.29592 51553 
1.28224 66152 
0.78105 14288 
0.78112 94887 
0.78794 34531 
0.79805 14791 
0.80914 20410 
0.84083 60622 
1.18672 06418 
1.16249 32414 
1.15341 42400 
1.14083 92066 
0.87774 86686 
1.10349 12379 
1.07216 02642 
1.06189 60132 
1.01435 60076 
0.99294 66336 
0.99575 27146 
1.00408 07498 


L. J. Comrie (1893-1950) 


Andoyer 


09199 
25798 
12322 
72317 
68504 
33564 
81363 
31849 
83411 
73164 
90501 
26311 
88324 
61502 
23502 
17498 
06214 
71500 
79501 
03055 
60463 
12500 
85730 
64518 
54770 
81500 
70743 
-57948 54234 
-59146 37910 
79624 
45886 
66590 
70350 


14178 
.63993 66167 
26347 
70538 


66536 
99574 27147 
07499 
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UNPUBLISHED MATHEMATICAL TABLES 
In this issue there is a reference to an unpublished table in RMT 1015 


147[A, P].—Lro Storcn, Admittance-Impedance Conversion Tables, Tech- 
nical Memorandum No. 274, Hughes Aircraft Co., Research and Devel- 
opment Laboratories, Culver City, California. 10 p. manuscript tabu- 
lated from punched cards. Copy deposited in UMT Fie. 


The table gives 4S values of (1 + s*)—! and s(1 + s*)— for s = 0(.001)1. 
It is intended to facilitate the calculation of the reciprocal of a complex 
number. The table is an extension of a table of JAHNKE & Empe. [4th ed., 
appendix, p. 13.] 


148[F].—F. GRUENBERGER, Lists of primes. Six sheets tabulated from 
punched cards. Deposited in the UMT Fie. 

These lists of primes are for the ranges 10100009 to 10132211 and 
50000017 to 50040013. The lists were computed on a CPC as a fill-in project, 
without attempting to program for speed. A graph showing the distribution 
of differences between consecutive primes in the ranges 1000003—1024523, 
10100009-10132211 and 50000017—50012839 is included. 


F. GRUENBERGER 
University of Wisconsin 
Madison, Wisconsin 


149[F].—A. GLopEN, Table of solutions of the congruence x +- 1 = 0 (mod p) 
for p < 20000. Manuscript, 2 p., deposited in the UMT Fine. 

P The table gives for each of the 16 primes p of the form 256k + 1 less 

than 20000, the 64 solutions of the congruence mentioned in the title which 

are less than $). 


A. GLODEN 


11 rue nee Jaurés 
Luxembourg 
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TECHNICAL DEVELOPMENTS 


AN AUTOMATIC COMPUTER IN AUSTRALIA 


An automatic computer, the “C.S.I.R.0. Mk. I Digital Computer,” 
designed and constructed by the Radiophysics Division of the Common- 
wealth Scientific and Industrial Research Organization in Sydney, Australia, 
is now in service. It is of the all-electronic, serial-binary type with a main 
store consisting of a group of ultrasonic delay lines with a total capacity 
of 1,024 words of 20 binary digits each. An auxiliary store in the form of 
an unsynchronized magnetic drum is incorporated with a capacity of 1,024 
similar words, later to be extended to 4,096 words. 

No attempt has been made to obtain very high speeds of operation, the 
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objective being to construct a computer which is simple and sufficiently 
flexible from the engineers’ and mathematicians’ point of view. Thus, 
equipment can be removed or added without dislocating the mode of 
operation. The device is intended at this stage primarily for research into 
numerical methods and programming techniques. 

Standard electronic components are used throughout the 1,500 or so 
tubes and circuits. The clock pulse rate is 333 kcs. A minor cycle of one 
word length covers a period of 60 usec, and the major cycle of 16 words, 
equal to the delay of each ultrasonic storage channel, covers a period of 1 
millisecond. On the average the time occupied by the selection of a com- 
mand and the corresponding operation is about 2 major cycles. 

Numbers and Commands. Numbers are registered in the straight binary 
scale, with each of 19 digits together with a sign digit in the most significant 
place, and with negative numbers being stored in complementary form. 
The convention with regard to the position of the binary point is only 
significant in the operation of the automatic multiplier, where it is placed 
immediately to the right of the sign digit of the total product. Although 
each number consists of only 20 binary digits, it is not difficult to programme 
for 40 digits or even greater accuracy. 

The scheme adopted for commands is of the ‘‘two address’ type, in 
which each operation is considered as a transfer of the content of one register 
to another, the former being specified as a “source” and the latter as a 
“destination.”’ An arithmetical function is specified as a quality of the par- 
ticular transfer demanded by a command. A command is divided into three 
groups of digits. Two adjacent groups of 5 digits define the source and 
destination, and a further group of 10 digits specifies a numerical component 
normally used to indicate which store position is called for when either of 
the other two addresses involve the store. If the store is not called, the third 
address may be used to store special information. 

These digit groupings are in the order in which digits are transmitted 
from a register: destination, 1—5; source, 6-10; numerical, 11-20. Commands 
are held in and are normally accepted from the store in serial order. 

Organization. The computer is of the serial type, and all transfers take 
place along a single ‘digit trunk.’’ All registers are connected to this trunk 
via ‘‘function gates” which are under the control of the central sequencing 
unit and the command decoding devices. The digit trunk consists of two 
parts, an “output trunk’ and an “input trunk,” and transfers are made 
between these conductors during single minor cycle periods, determined by 
a time selector which operates upon the detection of equality between digits 
11-14 of the current command and the current minor-cycle number. 

The selection of a single command from the store, and its performance, 
requires four transfers. The registers involved in controlling these actions 
are: 


(Ai) The “sequence register’: a 10-digit register which keeps a tally of 
the progress of the programme, and instructs whence the next command 
is to be withdrawn. Its contents are normally increased by unity following 
the selection of a command. 

(A ii) The “store control register,” of 10-digit capacity, which is con- 
nected in the upper 6-digit positions, i.e. (15-20), to a decoding selector 
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used to specify which store channel is needed, and in its lower 4 digit 
positions (11-14) to the time selector. 

(A iii) The “interpreter register,” of 20-digit capacity, which is con- 
nected to the ‘‘source and destination selectors” each possessing 32 outputs. 

Some of these registers possess special function gates and may be called 
to transfer or receive by a command, as follows: 

(B i) The sequence register can be called to “count in,” i.e., to increase 
its content by a number equal to the number of unit digits entering; to 
“add in” on the digit positions 11-20; and to “substitute in” on the same 
positions, i.e., reset and add in. Also it can be called to transmit in the digit 
group 11-20. 

(B ii) The interpreter register can be called to transfer out its numerical 
content on digit positions 11-20. 

The four transfers involved in satisfying a single command are of invari- 
able form. They constitute what is called the computer routine and are as 
follows: 

(Ci) The content of the sequence register is transferred to the store 
control register whose content it replaces. 

(C ii) Under the control of the content of the store control register, 
the store is allowed to transmit. The time selector finally allows the desired 
command to enter the interpreter register via the “input trunk.” 

(C iii) The command may call: for the store to transmit or receive, 
hence the numerical part of the command is transferred from the interpreter 
to the store control register, replacing any previous content. 

(Civ) Under the control of the time selector and the action of the 
source and destination selectors, a selected pair of function gates is actuated, 
and the desired transfer is performed. 


The function gates and their actions are listed in Table I. 

The Arithmetical Unit. The arithmetical unit consists of a group of 5 
ultrasonic delay-line registers of which A, B, C are of 20-digit capacity, 
whilst H is of 10 digits only. The fifth register D can store sixteen 20-digit 
numbers. Registers A, C, and D can add, subtract, and replace, whilst H 
can read in and replace in either of the digit groups 1-10 or 11—20. Register 
B is a non-adding register. 

Register A is the main accumulator and possesses a number of special 
gates for performing the functions of shifting right and left by one place, 
as well as reading out the most significant digit and the lowest digit. It is 
also capable of certain logical functions. Registers C and D are capable of 
reading out the sign digit of their content. This applies to any single number 
in D. Registers A, C, and D can read out their contents, and also may do 
so with simultaneous reset. 

Special functions associated with multiplication have been introduced 
which also may be used for other purposes. For automatic multiplication, 
registers A, B, and C are coupled together. Register C holds the multipli- 
cand whilst B receives the multiplier. The registers A and B are connected 
into series circulation together with an extra digit period, for a period of 
41 minor cycles during which the multiplier digits are investigated succes- 
sively and removed from circulation. The total product is built up by 
successive addition and shifting, until actually it occupies registers A and 
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TaBLE I. Function gates 
Command addresses: n; S, D 


Sources (S) Destinations (D) 
Read Out: Read In: 
High speed store position ‘‘n”’ 1) High speed store position ‘‘n”’ 
2) low speed store position ‘“‘n” 2) low speed store position ‘‘n”’ 
3) current command positions 11-20 3) shift to binary/decimal reading 
4) read current card column and shift (4) to output “print” 
to next column 
{} position of next command (5) to output and “punch” 
6) content of register A (6) sequence register and “‘count”’ 
7) content of A and clear tas sequence register and “substitute” 
8) content of A X 2 8) sequence register and ‘‘add” 
9) content of A X 4 (9) and reset into A 
10) sign digit of content of A (10) and add into A 
11) lowest digit of content of A 133 and subtract into A 
12) lowest unit of content of A ~ 0 12) and substitute A by conjunction of 
content and entry 
(13) content of B (13) and substitute A by disjunction of 
content and entry 
(14) content of B shifted one place to (14) to shift total product left one posi- 
right tion if unit digit received 
(15) sign digit of content of B shifted (15) and store in B 
one place to left 
16) content of C (16) and multiply in B 
17) sign digit of content of C (17) and reset into C 
18) content of D position ‘“‘n” +18} and add into C 
19) sign digit of content D position “‘n”’ 19) and subtract into C 
20) content of H on positions 1-10 (20) and reset into D position ‘‘n” 
21) content of H on positions 11-20 (21) and add into D position ‘‘n”’ 
22) one lowest digit (22) and subtract into D position ‘“‘n” 
23) content of hand-set register num- (23) H and store group 1-10 
(24) os 5 of hand-set register num- (24) H and store group 11-20 
r 


(25) and stop if unit digits received 


B, the multiplier being lost. It is arranged that the binary point of the 
product lies between the two most significant digits of the more significant 
component which lies in A. Corrections are automatically made to the 
product at the outset if either or both the multiplier and multiplicand are 
negative. 

For storage of double length products, B can read out one digit early, 
leaving its most significant digit free for the insertion of a sign digit if 
needed. Rounding off to 20 digits is performed by reading out the most 
significant digit of B with a delay of one digit period, i.e., as a digit in the 
last significant place, and adding it to A. Before rounding off it may be 
desired to retain as many digits as possible in one 20-digit word, so a “‘left- 
shift’’ function is available. This shifts the entire product one place to the 
left and may be called by a digit impulse supplied to its function gate. No 
right shift is available. 

The H register can also read out in either of the digit groups 1-10 or 
11-20 and so may be used for shifting groups of 10 digits to left or right by 
10 positions. 

Two special “constant registers” are provided. These are hand set to 
any desired set of digits for a calculation. They are particularly useful in 
maintenance tests and in testing programmes. 

A special constant gives a single digit output in the least significant 
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place (i.e., position 1), and a further destination is provided to stop the 
sequence. This function is performed if a non-zero digit passes this gate. 
A list of the sources and destinations is given in Table I. 

Use of the Functions. Whenever a command does not involve the store 
or register D, the digits in positions 11-20 may be used for other purposes, 
and special 20-digit constants may be compiled by transfer via the H 
register. Frequent use is found for special digit groupings in the positions 
11-20. 

Absolute transfer of control can be made by using the read-out faculty 
of the interpreter register and the read-in property of the sequence register, 
whilst relative transfers can be made by adding into the latter register. 
Conditional transfers are made by using the sign-selecting property of the 
arithmetical registers. A sign digit is counted into the sequence register, 
and an absolute or relative transfer may follow or not according to the sign. 
Multiconditional transfers can be made in a single count. 

For instance, in terms of the code numbers of Table I, the commands 


oe i @ 
a: 3, 3, 


respectively, replace the content of the sequence register by “mn” and add 
to it; whilst 
0; 10, 6 


adds a unit to the sequence register if the content of register A is negative. 

The use of standard sub-routines as aids to building complete pro- 
grammes is adopted. Most sub-routines can, with the aid of the properties 
of the sequence register and the multiplicity of adding registers in the 
arithmetical unit, be made independent of their actual position in the store. 
This simplifies programming. Transfers to and from sub-routines to the 
master programme are performed simply by use of the last word position 
of the D register, numbered 15. The content of the sequence register is 
stored there immediately before the transfer to a sub-routine headed, say, 
at position m. One of the commands of the sub-routine adds 1 to this address, 
and the last command transfers the content of Dy to the sequence register, 
which in turn calls for a transfer of control to the next command in the 
master programme. This method has the advantage of not involving the 
main arithmetical registers A, B, and C during transfers. The commands 
needed are shown in Table II. 

Input and Output. At present data are inserted via punched cards of the 
IBM type. These are read in a columnar fashion, ten binary digits being 


Tas_e II 
Command Position Command Significance 
n—2 15; 5, 20 Stores 2 — 1 in Dy 
n—1 m; 3,7 Transfers to “‘m”’ 
m _ Head of sub-routine 
m+r-—1 15; 22, 21 Adds unit to Djs 
m+r 15; 18,7 Last command of sub-routine; 


transfers to “n” 
n _ Continues programme 
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entered from each column. A primary routine of twenty commands is 
entered, prior to the feeding of cards, through a group of permanently wired 
stepping switches. This routine supplies sufficient information to the com- 
puter to enable it to assemble 10-digit data into 20-digit words. 

The card reader is also capable of reading decimal data, punched in the 
usual decimal manner, and changes from binary to decimal reading and 
back are by the computer. 

Most programme data are supplied in binary form being mostly of sub- 
routine type, and little extra effort is required to punch the cards of the 
master programme in binary form. 

The X and Y punch positions are used for reading control; in particular 
a Y punch is used to adjust commands as they enter the computer in order 
to allow for data which vary with the problem and so on. 

Output is obtained in two ways: by standard teletype page printing 
and by decimal card punching, again in columnar fashion. Punching is so 
arranged that results may be re-inserted later into the computer as problem 
data if required. 

T. PEARCEY 
Division of Radiophysics 
Commonwealth Scientific and Industrial 
Research Organization 
Australia 
DISCUSSIONS 


Minimum Access Programming 


The following discussion outlines programming techniques which may 
be applied to reduce computation time on certain large scale general purpose 
computers. The computers to which these techniques can be applied most 
readily have memory systems of large capacity, storing commands and 
operands in an alterable medium. These techniques have been applied to great 
advantage in the ERA 1101, a computer using a magnetic drum storage 
system. On the 1101 the techniques have given an advantage of 20 to 1 in 
computing speed over ordinary induction loop techniques for the average 
program, with an advantage of 100 to 1 over short stretches. 

The internal memory of a typical computer under discussion is divided 
into compartments called memory boxes, each memory box having an 
address. The computer performs operations such as addition and subtraction 
on operands contained in memory boxes. The content of a memory box is 
called a ‘‘word.”” The word may either be an operand or a command.! If 
the word is a command, it will consist of a command code and one to four 
execution addresses. 

The addresses of the memory boxes have a cyclic order. Access to a 
memory box is made possible when a coincidence is detected between a 
locating register and a storage address register. The storage address register 
contains the address of the memory box to which access is desired. As each 
memory box is scanned its address is held in the locating register. The 
addresses of all the memory boxes pass through the locating register in 
their cyclic order every memory cycle. 

Let an address be sent to the storage address register. The interval 
between the time that the coincidence detector is turned on and the time 
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that a coincidence is detected depends on the distance between the memory 
box containing the instruction and the memory box containing the operand. 

For each command of a computer there is a minimum interval between 
the location of the command and its operand or operands, such that if a 
different interval is programmed the effect would be a loss of operating 
time. Attention to these minimum (coding) intervals in construction of a 
program is called minimum access programming. 

In order to examine minimum access techniques more closely, those 
used on a specific machine, the ERA 1101, will be discussed. 

The ERA 1101 is a high-speed electronic, general purpose digital com- 
puter, having a large internal memory. The computer operates on 24 binary 
digits (bits) and uses the one address system of logic. The internal memory 
consists of a magnetic drum? of 16,384 memory boxes. A portion of the 
drum in which a digit is stored is called a memory cell and a portion of 
the drum in which a word is stored is called a memory box. For each memory 
box there are 24 memory cells. All of the digits of a word are read from or 
written in a memory box in parallel. Arithmetic operations treat numbers 
up to +2 — 1 in a double length accumulator. 

The drum is subdivided into eight groups along its length, so that for 
each angular position of the drum, any one of eight memory boxes is avail- 
able for access. Each group, by: itself, contains 2,048 memory boxes. 

Let a memory box on the drum be designated by the address aj302 - - - 
20,4 held by locating registers, where the a’s are binary digits. Then the 
digits @,4yd) in the locating registers denote the group containing the 
memory box, the digits a:3¢:.2 denote the quadrant containing the memory 
box, and the remaining digits denote the angular position of the memory 
box within the quadrant. 

































































Fic. 1. Address mapping on drum for 1 interlace. 
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Fic. 2. Patching connections for 1B interlace. 


In order to select memory locations around the drum so that they are 
associated with consecutive memory boxes, the storage address register 
(SAR) is connected to the locating register (LR) through an interchange- 
able plugboard called an interlace chassis. 

If the digits of the SAR are carried into their corresponding positions 
in the LR, then consecutive addresses are located in consecutive memory 
boxes. This is called a one interlace. 

For adjacent memory boxes, with any interlace, the difference between 
their addresses is called the increment. The one interlace has an increment 
of one. 

One may see that there are 14! ways to interlace the drum with the 
storage address register, although only a few of these have been useful. 
Those used most frequently are the 1, 2, 4, 8, 16, 32, 64, 128 and 1B inter- 
lace. The special interlace hook-up called the 1B interlace has been found 
useful for loading the drum. The interlace may be dialed. 

With a 64 interlace, memory boxes having consecutive addresses are 64 
cell periods apart. One cell period is the interval between two adjacent 
memory boxes on the drum and is equal to 8 microseconds. Also for a 64 
interlace we have an increment of 32. 

Starting with address zero, adjacent memory boxes have the following 
addresses for 64 interlace. Octal notation is used. 


X is the group index: 


0X000 
0X040 
0X 100 
0X 140 





0X740 
1X000 
1X040 
1X100 





1X700 
1X740 
2X000 
2X040 
2X 100 





etc. 
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Each command consists of a six bit command code, dz, ---, dis, four 
skip bits a7, ---, @4 and a fourteen bit execution address a;, ---, do. 

A 1 appearing in any one of the skip bit positions of a command can be 
used to advance the program address counter by pulsing one of its first ten 
stages. Four selector switches determine which stages are to be pulsed by 
the skip bits should they appear. 

Delays inserted in the skip bit lines allow multiple skip bits to be used 
with the same command. 

One skip bit setting used by the author carries 


ds; into 2°, dys into 2', and ay into 2? 


with the selector switch for ay turned off. 

Then if a command appears in box y with 1’s in ais, a1. and ay, the next 
command will appear in address y + 2° + 2! + 2? + 1. 

If an interlace of one were used, with no skips, at least one drum revo- 
lution would be required for each instruction. By choosing the proper inter- 
lace and skips, however, operands and instructions may be so placed that 
several instructions are executed in one drum revolution. 

In order to coordinate the philosophy of minimum access programming 
with techniques used in the ERA 1101, a discussion of the 1101 commands 
is in order. 

The list given in Table I gives the commands and their corresponding 
minimum allowable access times. In the table interpret ‘‘y’’ as the address 
of the memory box containing the word ‘“(y).’’ The symbol “‘—”’ is inter- 
preted as “goes to.”’ AL is the left half or higher order accumulator and 
AR is the right half or lower order accumulator. A is a 48-bit electronic 
register and Q is a 24-bit electronic register. Negative numbers are expressed 
by complements. 

The intervals A, B, C and D specified in Table I are minimum values. 
Should a shorter interval be programmed, the only effect is a loss of oper- 
ating time. The interval is then automatically lengthened by the duration 
of a complete drum revolution. 

The intervals of commands 17, 27, 35 and 36 are valid for a single 
writing operation. The minimum allowable interval between two successive 
writing operations in the same group is 256 cell periods. In different groups, 
the minimum interval is 4 cell periods. Note that a reading operation will 
normally occur between successive writing operations. 

The minimum interval for either C or D of commands 25 and 26 is 
approximately 6 + 0.33% cell periods where k is the number of places 
shifted. 


2'5 2'2 oll 2! 2? 2° 2” 2° 2° “ 2 22 2 2° 


Fic. 3. Patching connections for a 64 interlace. 
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TaBLeE I. List of commands and minimum time intervals for coding 


Address of Address of 
Address of Operand for Address of Operand for 
Order 2 Order n Order n + 1 Order n + 1 




















Z 
~ 
D 
— 
Minimum Allowable 
Command Time in Interval 
Symbol Code Command in Order in Cell Periods 
A B Cc. Ps 
y—>cAt 11 Clear Add: Insert (y) in A. 4 4 —- — 
y—>hAt 12 Hold Add: Add (y) to (A). 4 4 —-— — 
y—>cA- 13 Clear Subtract: Insert the complement 4 4—->- — 
of (y) in A. 
y—> hA- 14 Hold Subtract: Add the complement 4 4—-> — 
of (y) to (A). 
y>Q 16 Load Q: Insert (y) in Q. 4 4—-— — 
QaAmy 17 Substitute Digits: Insert the lower half 7 s—- — 
of (A) into y where there are ones in 
Q and block any change in y where 
there are zeros in Q. 
y—> mcAt 21 Absolute Clear Add: Insert the abso- 4 4- => 
lute value of (y) to A. 
y— mhAt* 22 Absolute Hold Add: Add the absolute 4 4—- — 
value of (y) to (A). 
y—> mcA- 23 Absolute Clear Subtract: Insert the 4 4o-_ — 
complement of the absolute value of 
(y) in A. 
y—> mhA- 24 Absolute Hold Subtract: Add the com- 4 4-—-> — 
plement of the absolute value of (y) 
to (A). 
ALk 25 Shift A Left: Shift A circularly, the _-_ - - = 
number of oon designated by the 
number in the execution address part 
of the instruction. 
QLk 26 Shift Q Left: Shift Q circularly, the -_ - - = 
number of — designated by the 
number in the execution address part 
of the instruction. 
Ap—>y 27 Substitute Execution Address: Insert 4 s;s—- — 
(AR) into y blocking any change in 
the upper six binary digits of (. 
y—scAt 31 Split Clear Add: Insert (y) in A mul- 4 4 —- — 
tiple precision. 
y— shAt 32 Split Hold Add: Add (y) to A multiple 4 4 —- — 
precision. 
y—> scA- 33 Split Clear Subtract: Insert the com- + 4 —>. — 
plement of (y) in A multiple precision. 
y—> shA~ 34 Split Hold Subtract: Add the comple- 4 4—-> — 


ment of (y) to (A) multiple precision. 
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TaBLeE I—Continued 


Minimum Allowable 
Command Time in Interval 
Symbol Code Command in Order in Cell Periods 
A-y 35 Store A: Insert (AR) into y. 4 s—_- — 
Q-y 36 Store Q: Insert (Q) into y. 4 s—- — 
zero—> AR 37 Clear AR: Clear the lower order accu- —- —- 5 — 
mulator. 
Q—cA 41 Clear Add from Q: Insert (Q) in A. —- —- s&s — 
Q—hA 42 Hold Add from Q: Add (Q) to (A). —- —- s&s — 
A—~cQ 43 Transmit A to Q: Insert (AR) in Q. —_ —- s— 
Qs 44 Q Jump: If (Q) is negative, insert (y) —- —- s — 
in the SAR. In either case shift Q 
circularly one place to the left. 
J 45 Jump: Insert y in the SAR. —_ —- 5s — 
cJ 46 Sign-Conditional Jump: If (A) is nega- —- —- 5s — 
tive, insert y in the SAR. 
zero J 47 Zero-Conditional Jump: If (A) is zero — — 5§ — 


insert y in the SAR. In either case 
subtract one from A. 


Qy—>cA 51 Clear Logical Multiply: Insert (y) in 4 4—- — 
AR where there are ones in Q, block- 
ing a change in AR where there are 
zeros in Q. 


Qy—hA 52 Hold Logical Multiply: Add (y) to 4 4—> — 
(AR) where there are ones in Q, block- 
ing a change in AR where there are 


zeros in Q. 

y>P 53 Print Only: Print the six bit character + 4—> — 
whose address is y. 

y—>P/P 54 Print and Punch: Print and punch the 4 4 —-_ — 
six bit character whose address is y. 

Is 55 Intermediate Stop: Stop the computer - - -—- =— 


and turn on the intermediate stop 
light. Proceed only after the operator 
pushes the start button. 


Os 56 Optional Stop: Stop the computer if -_ - =—- = 
the optional stop button has been 
pushed, and give suitable indication. 


FS 57 Final Stop: Stop the computation and -_ - -—- =— 
turn on the final stop light. 
Q-y—>cA 61 Clear Multiply: Obtain (y)-(Q) in A. 4 4— 37 
Q-y—>hA 62 Hold Multiply: Add (y)-(Q) to (A). 4 13 — 4§ 
A/y—'cQ 63 Divide: Obtain (A)/(y). Find the quo- 4 4— 4 
tient in Q and the non-negative re- 
mainder in A. 
65 Optional J ~ Same as optional a —-_ —- 5s — 
except when the computer ison TEST. 


If on TEST, and the optional jump 
switch is on, then this command is 
the same as the jump command. 


66 Do Nothing: Performs no operation. —- — s — 
yt+1—~A 71 par AOS ee ae Insert (y)+1 in 4 4 —>— — 











178 AUTOMATIC COMPUTING MACHINERY 


The intervals of commands 53 and 54 are for a single print operation. 
The minimum allowable interval between two successive print operations 
is approximately 19,000 cell periods or 150,000 microseconds. Should a 
shorter interval be programmed it will automatically be lengthened to this 
interval. Note, however, that while a print command is being executed, the 
computer can proceed with the program of instructions provided the follow- 
ing orders do not call for another print operation. 

Although interval B of command 62 is indicated as being a minimum of 
13 cell periods, a 12 cell period may be used for B when interval A is pro- 
grammed at no more than 4 cell periods. Under this condition, interval C 
is therefore 16 cell periods. 

The 66 command performs no operation but its skip bits may be used 
to augment those of the previous command. 

In writing a minimum access program, the programmer sometimes takes 
advantage of induction loop programming methods. A case in point is the 
matrix multiplication of two 10 X 10 matrices. The solution of the problem 
gives X - Y = Z, where X, Y and Z are 10 X 10 matrices. 

A portion of the program of a matrix multiplication is shown in Table II. 
Skips are not used as the program was written before they were made 
available. In preparation for the execution of this program, elements of X 
are inserted row by row in memory boxes whose addresses are designated 
by the following table of addresses at 128 interlace in group ¢. Octal nota- 
tion is used. 


401 402 403 404 405 406 407 410 411 412 
421 422 423 424 425 426 427 430 431 432 
441 442 443 444 445 446 447 450 451 452 


621 622 623 624 625 626 627 4630 631 + ~# 632 
Table of Addresses 


The elements of Y transpose are inserted row by row in memory boxes 
whose addresses are designated by the table of addresses at 128 interlace, 
but in group 1. 

Z is written in a similar block, where the elements of the first row are 
in consecutive addresses at 64 interlace. 

The program of Table II obtains the products of the rows of X and a 
column of Y in 30 revolutions of the drum, and the entire matrix product 
has been clocked at 8.9 seconds where no overflow occurred during the 
multiplication. If a 2;; occurs larger than 2* — 1, then 2% — 1 is subtracted 
from 2;; and a conditional jump is executed. 

The scaling routine which starts at address 30540 does not utilize mini- 
mum access techniques. This is an ordinary induction loop routine. A 
minimum access attack may be made on the scaling routine at the expense 
of memory space. 

It is possible to manufacture portions of the minimum access program 
in the computer. This manufactured program may be punched out ready 
for use. 

From the program of Table II and from other material presented above, 
it may be concluded that: Arranging operands and instructions so that they 
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TABLE II. Matrix multiplication on 64 interlace of two 10 X 10 matrices 
Insert the following commands at 1 B interlace (Octal notation is used): 


Memory Box 


03120 
03121 
03122 
03123 


03131 


03220 
03221 


03231 


03320 
03321 


03330 
03331 


At 64 interlace: 
30500 


30501 
30502 
30503 
30504 
30505 
30506 
30507 


30534 
30535 
30536 
30537 


31001 
31002 
31003 
31004 


31023 
31024 
31025 
31031 


31037 


i1 
27 
27 
27 
45 
45 
45 
45 


Command 
03 10 
03 «10 
03 «(11 
03° «11 
03 «14 
03 «(00 
03 «00 
03 «(04 
01 00 
01 01 
01 04 
03 «05 
03 05 
03 «05 
03 «05 
03 «05 
03 05 
03 «(05 
03 «05 
00 00 
00 00 
03 07 
77 «77 
00 00 
01 10 
01 30 
01 10 
01 40 
01 10 
01 70 
03 «00 
01 20 
00 00 


02 


Explanation 


3120 at 1B interlace be- 
comes 10031 at 64 interlace. 


The operation starts with 
the command in memory 
box 10031. 


3220 at 1B interlace be- 
comes 10032 at 64 interlace. 


3220 at 1B interlace be- 
comes 10033 at 64 interlace. 


Minus ten counter. 
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Subroutines starting in memory boxes 31041, 31101, 31141, 31201, 31241, 31301, 31341, 
and 31401 have been omitted to conserve space. Their function is similar to the sub- 


routines of 31001 to 31025 and 31441 to 31465. 


31441 
31442 
31443 
31444 


31464 
31465 


31471 





14 
30 
14 
40 
70 
00 
24 


42 
03 
44 
05 
25 
00 
40 
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TABLE II1—Continued 


Memory Box Command Explanation 
30000 1i O01 00 O1 zero A 
30001 16 01 00 02 wu Q 
30002 62 O01 30 03 QO-nj—>hA 
30003 16 O01 00 04 X12—> QO 
30004 62 01 40 05 Q-923-—> hA 
30023 16 O01 00 24 X10 QO 
30024 62 O01 70 25 Q:9¥10;3 > hA 
30025 46 03 20 26 
30026 14 O01 20 27 37777777 — hA- 
30027 46 03 00 31 
30030 45 03 05 40 
30031 12 O1 20 32 37777777 — hAt 
30032 35 00 00 00 
30033 45 03 00 40 


The subroutines starting in memory boxes 30040, 30100, 30140, 30200, 30240, 30300, 
30340, and 30400 are similar to the subroutines of 30000 to 30033 and 30440 to 30472. 


30440 11 O1 00 O1 zero —>cAt 
30441 16 O01 04 42 %10,1-> Q 
30442 62 01 30 03 O-n3—>hA 
30443 16 O1 04 44 x0,2—> Q 
30444 62 O01 40 05 Q-23—> hA 
30463 16 O01 04 64 0,0 Q 
30464 62 O01 70 25 Q-910; > hA 
30465 46 03 24 66 
30466 14 01 20 27 37777177 > hA 
30467 46 03 04 71 
30470 45 03 05 40 
30471 12 01 20 32 
30472 35 00 00 00 
30473 71 03 05 37 
30474 35 03 O05 37 
30475 46 01 30 30 
30476 57 00 00 00 

At 1 B interlace: 
12660 12 01 20 27 12660 at 1 B interlace is 
12661 12 01 20 27 32026 at 64 interlace. 
12671 12 O01 20 27 
12760 46 03 06 00 12760 at 1 B interlace is 
12761 46 03 06 00 32027 at 64 interlace. 
12771 46 03 06 00 
13060 14 O1 20 £31 13060 at 1 B interlace is 
13061 14 01 20 31 32030 at 64 interlace. 
13071 14 01 20 31 
30600 25 00 00 31 
30601 43 00 00 00 
30602 11 03 OS 75 
30603 35 03 OS 74 


30604 44 03 06 06 








Memory Box 


30605 
30606 
30607 
30610 


30636 
30637 


At 1 B interlace: 


17120 
17121 


17131 


17220 

iia 
17222 

00, 

72. 17223 


17231 


17320 
17321 
17322 
17323 


17330 
17331 


At 64 interlace: 


31731 
31771 


13030 


13160 
13161 
13162 
13163 


13170 
13171 


At 64 interlace: 


12027 
12030 
12031 
12032 


30540 
30541 
30542 








Command 
03 «05 
03 05 
03 «(05 
03 06 
01 10 
00 00 
03 17 
03 17 
03 17 
03 00 
03 00 
03 «(01 
03 «OO 
03 04 
01 30 
01 31 
01 31 
01 32 
01 34 
03 «(05 
00 00 
00 00 
03 00 
03 «00 
03 00 
03 «Ol 
03 «(O1 
03 04 
03 «(04 
77 «(17 
7 OF 
7 «697 
77 «177 
00 00 
00 00 
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TaBLe II1—Continued 


50 
74 
74 
04 


02 
30 


31 
71 





Explanation 


One — hAt* 
40 —» hAt 


17120 at 1 B interlace is 
13031 at 64 interlace. 


17220 at 1 B interlace is 
13032 at 64 interlace. 


17320 at 1 B interlace is 
13033 at 64 interlace. 


13160 at 1 B interlace is 
32031 at 64 interlace. 


This is the scaling routine. 
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TaBLe II1—Continued 


Memory Box Command Explanation 
30554 35 03 OS 73 
30555 11 03 06 76 
30556 27 03 OS 60 
30557 27 03 OS 62 
30560 11 00 00 00 
30561 25 00 00 00 
30562 35 00 00 00 
30563 11 03 OS 57 
30564 so S&S Bw HH 
30565 27 03 OS 60 
30566 27 03 OS 62 
30567 14 03 O05 72 
30570 46 03 05 60 
30571 455 01 00 31 
30572 16 O1 14 64 
30573 00 00 00 00 
30574 00 00 00 00 
30575 "7 um hu 
30576 77 «77 7%6©«©33)0)=«6~——ten? 
30577 00 00 00 00 


occur as close as possible to the minimum allowable coding intervals shortens 
the computation time considerably. In the program of Table II the method 
used to obtain a product of a row by a column is faster than an ordinary 
induction loop program by a factor of approximately 50. A 10 X 10 matrix 
multiplication may be performed in approximately 8 seconds using minimum 
access techniques. 


1In the UNIVAC a word contains two commands. 
2 A. A. CoHEN, ‘Magnetic drum storage for digital information processing systems,” 
MTAC, v. 4, p. 31-39. 
Engineering Research Associates 
Arlington, Virginia 
Davin P. PERRY 
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This is a comprehensive exposition of the logical organization of auto- 
matic digital computers and of the principles of their physical realization. 
It is the first publication of its kind in German and far better than most 
similar writings in English. Concisely and precisely written, complete, and 
as unbiased as is possible in this controversial field, it is highly recommend- 
able as an introduction to the subject. 

The first chapter deals briefly with the why and how of automatic 
computation, its history, relations to formal logic and neurophysiology, 
and (perhaps too briefly) applications. A second chapter lists and describes, 
in only four pages, the major components of automatic computers. The 
third chapter explains in full detail the various methods for representing 
numbers (binary, coded decimal; fixed vs. floating point) and for carrying 
out the elementary arithmetic operations. There follows a chapter on coding 
and programming, including such topics as flow diagrams, branch points, 
single and multiple address codes, modification of instructions by arith- 
metic operations, and a completely carried out example (multiplication of 
two matrices). Methods of evaluating mathematical functions as well as 
checks, both built-in and programmed, are also discussed. The final chapter 
describes the basic circuits occurring in computers, such as ‘‘and’’ gates, 
“or” gates, adders, etc.; the several methods of number storage (acoustic, 
magnetic, electrostatic); and media for input and output. There is a sum- 
mary table of computer developments as of December 1949 and a compre- 
hensive bibliography. 

FRAnz L. ALT 

NBSCL 

NEws 


Association for Computing Machinery.—The spring meeting of the Association was 
held on May 2 and 3, 1952, at the Mellon Institute, Pittsburgh, Pennsylvania. On May 2nd 
at 7:00 p.m. a banquet was held at which time a talk, ‘History of mechanical computing 
machinery,” was presented by GEorGE C. CuasE, Monroe Calculating Machine Company, 
Orange, N. J. The program for the meeting was as follows: 


May 2, 1952, 10:00 a.m. 


General Session Franz L. Att, NBSCL and President of 
ACM, Chairman 

Address of Welcome E. R. WEIDLEIN, President, Mellon Insti- 
tute 

Evaluation of Automatic Computing IrvEN Travis, Burroughs Adding Machine 


Co. 





arte 
102 


ito- 
on. 
ost 
and 
nd- 


itic 
gy, 
es, 
The 
ing 
ing 
ing 
its, 
ith- 
. of 


yter 
tes, 
tic, 


re- 


was 
2nd 
ting 
ny, 


t of 
sti- 


hine 





AUTOMATIC COMPUTING MACHINERY 187 


Transistor Physics 
Mark IV 


W. SHocKLEY, Bell Telephone Laboratories 
H. H. Arxen, Harvard Computation Lab- 
oratory 


May 2, 1952, 2:30 to 5:00 p.m. 


Computers and Components I 

Special-Purpose Digital Data-Processing 
Computers 

The Remington-Rand Punch-Card Elec- 
tronic Computer, Type 409-2 

The Elecom 100 General Purpose Computer 


The Quadratic Arc Computer (QVAC) 

A System for Counting and Recording Elec- 
trical Impulses in Printed Decimal Form 

The Logical Organization of the New IBM 
Scientific Calculator 


Jan. A. Raycuman, RCA. Chairman 

B. M. Gorpon and R. N. Nicora, Labora- 
tory for Electronics, Inc. 

L. P. Crosman, Remington-Rand Corp. 


ALBERT AUERBACH, Electronic Computer 
Corporation 

M. J. MENDELSON, Northrup Aircraft, Inc. 

J. L. Linpesmita, Clary Multiplier Corp. 


N. Rocuester, IBM Corp. 


May 2, 1952, 2:30 to 5:00 p.m. 


Mathematical Applications I 

Engineering Problems Requiring Automatic 
Computing Facilities 

Digital Computer Methods for Problems 
which Involve Linear Inequalities 

Computational Problems of Linear Pro- 
gramming 

Small Problems on a Large Computer 

Firing Table Computations on the ENIAC 

Small-Scale Research and Automatic Com- 
puting 


Mina REEs, ONR, Chairman 
E. L. Harper, Westinghouse Electric Corp. 


ALEX Orpen, Ha. U.S.A-F. 


A. Carnes, Carnegie Institute of Tech- 
nology 

C. W. Apams, MIT 

H. L. REED, Jr., Aberdeen Proving Ground 

E. C. BERKELEY, Edmund C. Berkeley & 
Associates 


May 2, 1952, 2:30 to 5:00 p.m. 


Theory I—Information and Control 

Use of Computing Machinery in Applica- 
tions of Information Theory 

An Upper Bound on the Informational Ca- 
pacity of a Synapse 


Automatic Control of Machinery 


The Maze Solving Computer 


A Chart Method for Simplifying Truth 
Functions 


J. H. Curtiss, NBS, Chairman 
W. G. TuLLer, Melpar, Inc. 


W. S. McCu..ouca, Neuropsychiatric In- 
stitute, DonaLD MacKay, Kings College, 
London 

P. L. Nres, G. P. Tanguary, D. R. Aur- 
DERHIEDE, D. W. Brown, C. J. Jacosy, 
R. B. Ketter, Harvard Univ. 

R. A. Watvace, L. and O. Research and 
Development Co. 

E. W. Vertcu, Burroughs Adding Machine 
Co. 


May 3, 1952, 9:30 to 12:00 Noon 


Computers and Components II 


Transistor Circuits for Computers 

Standardized Printed Circuit Units for Digi- 
tal Computers 

Non-Linear Switching Elements 


Joun D. Ditton, Air Force Missile Test 
Center, Cocoa, Florida, Chairman 

J. H. Fevker, Bell Telephone Laboratories 

D. L. Jonnson, Elliott Brothers, Ltd., Bore- 
hamwood, England 

B. Morrat, F. A. ScHwertz, Mellon Insti- 
tute, B. O. MarsHa.L, Air Force Cam- 
bridge Research Center 
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Optical Elements for Computers 


The Selenium Rectifier—A Non-Linear and 
Asymmetric Resistance Element 


B. O. MarsHati, Air Force Cambridge 
Research Center, J. R. Bowman, F. A. 
ScHWERTZ, Mellon Institute 

N. Harpy, International Resistance Co. 


May 3, 1952, 9:30 to 12:00 Noon 


Mathematical Applications II 


Discussion on the Use and Construction of 
Subroutines 


Aston S. HousEHOLDER, Oak Ridge Na- 
tional Laboratory, Chairman 

JaMEs ALEXANDER, Argonne National Lab- 
oratory, HERMAN H. Goxtpsting, IAS, 
Josern H. Levin, NBS, H. RuBENSTEIN 
and J. D. RUTLEDGE, Remington-Rand, 
Inc. 


May 3, 1952, 9:30 to 12:00 Noon 


Theory II—Mathematical 

A Unified Approach to the Monte Carlo 
Method 

The Solution of Boundary Value Problems 
by the Method of Kernel Functions 

Boundary Value Problems in Doubly Con- 
nected Domains 


C. V. L. Smita, ONR, Chairman 
J. H. Curtiss, NBS 


STEFAN BERGMAN, Stanford Univ. 


Franz L. Att, NBS 


May 3, 1952, 2:00 to 5:00 p.m. 


Computers and Components III 


Digital Storage Using Ferromagnetic Ma- 
terials 


Some Recent Research on Ultrasonic Propa- 
gation in Solid Media 

Static Magnetic Memory—Its Application 
to Computers and Controlling Systems 

Static Magnetic Memory for the ENIAC 


Magnetic Binaries in the Logical Design of 
Information Handling Machines 


L. P. CrosMAn, Remington-Rand Corpora- 
tion, Chairman 

P. D. Atkinson, A. E DEBarr, R. MILLER- 
sip, R. C. Rosstns, Elliott Brothers, 
Ltd. ' 

T. F. Rocers, Air Force Cambridge Re- 
search Center 

An WANG, Wang Laboratories 


Isaac L. AUERBACH, Burroughs Adding 
Machine Co. 
N. B. SaunpERs, Transducer Corp. 


May 3, 1952, 2:00 to 5:00 p.m. 


Mathematical Applications IV 


Discussion on the Use and Construction of 
Subroutines 


Atston S. HousEHOLDER, Oak Ridge Na- 
tional Laboratory, Chairman 

RosELyn Lrpxis, NBS, Davip J. WHEELER, 
University of Illinois and University of 
Cambridge, Joun W. Carr, MIT, GRACE 
M. Hopper, Remington-Rand, Inc. 


May 3, 1952, 2:00 to 5:00 p.m. 


Theory III—Logic and Circuit Synthesis 


Formal Logic and Switching Circuits 


Minimization Synthesis of Two-Valued 
Feed-Back Circuits 


H. H. Arxen, Harvard Computation Lab- 
oratory, Chairman 

THEODORE Katin, Air Force Cambridge 
Research Center 


Wittiam Burkwart, Monroe Calculating 
Machine Co. 
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Characteristic Numbers and Their Use in 

the Decomposition of Switching Functions 
Rectifiers as Elements of Switching Circuits 
The Theory of Counting Techniques 


The Application of Counting Techniques 


Warren L. Semon, Harvard Computation 
Laboratory 

Peter F. Stronc, Harvard Computation 
Laboratory 

THEODORE SINGER, Harvard Computation 
Laboratory 

Rosert L. AsHENHURST, Harvard Compu- 
tation Laboratory 


University of California at Los Angeles.—On April 30, May 1-2, 1952, at UCLA, a 
symposium on electronic computers was held under the sponsorship of the Los Angeles 
IRE Professional Group on Electronic Computers and the Department of Engineering of 


the University. The program was as follows: 


Wednesday, April 30, 1952 


Registration 
Opening Session 
Introduction 


Welcome 


Keynote—Engineering Tomorrow's Com- 
puters 
Session on Magnetic Devices 


Design Features of a Magnetic Drum Mem- 
ory for the National Bureau of Standards 
Western Automatic Computer (SWAC) 

Problems Involved in Magnetic Tape Re- 
cording 

Survey of Tape Drive Systems 


Session on Analog Devices 


An Electro-Mechanical Multiplier for Ana- 
log Computer Application 

The Thermal Analyzer—A Special Purpose 
Analog Computer 

Panel Discussion—Utilization of Germa- 
nium Diodes 


Panel Discussion—Designing for Maximum 
Reliability 


R. G. CanninG, Naval Air Missile Test 
Center, Point Mugu, Calif. 

L. M. K. BoeE.tter, Department of Engi- 
neering, UCLA, Chairman 

H. D. Husxey, NBSINA 


Joun J. ConNnoLLy, The Ranp Corporation, 
Chairman 
R. THORENSEN, NBSINA 


Norman E. Grsss, Raytheon Manufactur- 
ing Company 

HAROLD SaRKISSIAN, Computer Research 
Corporation 

W. L. Martin, Department of Engi- 
neering, UCLA, Chairman 

S. E. Dorsey, Naval Ordnance Test Sta- 
tion, Inyokern, Calif. 

W. L. Martin and R. BromsBerc, Depart- 
ment of Engineering, UCLA 

L. L. Kitpatricx, North American Avia- 
tion, Inc., Chairman 

Norton BELL, Consolidated Engineering 
Corporation 

L. S. Petrrey, Hughes Aircraft Co. 

W. Speer, Computer Research Corp. 

A. S. Zuxin, Hughes Aircraft Co. 

Harry T. Larson, Hughes Air-craft 
Company, Chairman 

Joun J. Connoty, The Rand Corp. 

Harry D. Huskey, NBSINA 

Ropert Lusser, Naval Air Missile Test 
Center, Point Mugu, Calif. 

Rosert Raw.ins, Lockheed Aircraft Cor- 
poration 

Wittram WaAGENSEIL, Hughes Aircraft 
Company 
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Thursday, May 1, 1952 


Session on Programming and Coding 

An Approach to the Use of the IBM Card- 
Programmed Electronic Calculator in the 
Solution of Engineering Calculations 

Some General Precepts for Programmers 

Programming for On-Line Computations 

The Human Computer’s Dreams of the 
Future 

Automatic Program Control Utilizing a 
Variable Reference for Addressing 

Programming for Finding Characteristic 
Values of Mathieu’s Equation and the 
Spheroidal Wave Equation 

Session on Input-Output Equipment 


The Benson-Lehner Photoformer 

Input-Output on the New IBM Scientific 
Computer 

An Accurate Digital-Analog Function Gen- 
erator 

Some Techniques of Analog to Digital Con- 
version 

The Teleplotter, A Digital Plotting Device 

Summary Session 

Introduction 

Speakers 


RosELyn Lrpxis, NBSINA ,Chairman 
Murray L. LEsseEr, Northrop Aircraft, Inc. 


E. C. YowEtt, NBSINA 

H. LuxENnBERG, Hughes Aircraft Co. 

Ipa RuHopEs, NBS (Read by R. R. ReEy- 
NOLDs, NBSINA) 

A. S. Zuxtn, Hughes Aircraft Co. 


GERTRUDE BLancH, NBSINA 


H. Doeveman, Electronic Engineering 
Company of Calif., Chairman 

D. L. Prrman, Benson-Lehner Corp. 

M. M. Astranan, IBM Corporation 


W. A. FarrAnp, North American Aviation, 
Inc. 

H. E. Burke, Consolidated Engineering 
Corporation 

DonaLp F. BELLorr, Telecomputing Corp. 


R. L. Stsson, Computer Research Corp. 

D. H. LeEnmer, NBSINA 

Louis N. Ripenour, International Tele- 
meter Corp. 


Friday, May 2, 1952 
Tours were conducted to the following facilities: 


The Rand Corporation, 
Telecomputing Corporation, 
California Institute of Technology, 
Northrop Aircraft, Inc. 


OTHER AIDS TO COMPUTATION 
Linear Algebraic Systems and the REAC 


1. Introduction. A great variety of problems in both pure and applied 





mathematics involves, either directly or indirectly, the solution of systems 
of simultaneous linear algebraic equations. Although the solution of such a 
system can readily be indicated by determinants, it is found that the attain- 
ment of actual numerical answers frequently becomes laborious for systems 
of order greater than three. Further, it is found that a great deal of effort 
has been, and is being, expended in the development of numerical procedures 
and in the design and development of computers which can be applied to 
such linear systems. 

The approach in this instance is somewhat different. Here we are inter- 
ested in extending the utility of an existing computer, or to be more specific, 
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we are interested in the application of an electronic differential analyzer to 
linear algebraic systems. The methods herein presented should be applicable 
to various computing machines of this type. The specific computing machine 
referred to in the following pages is the REAC (Reeves Electronic Analogue 
Computer). 

The discussion immediately following this introductory section is a brief 
(and perhaps oversimplified) description of REAC techniques needed to 
instrument the methods proposed herein. Thereafter, certain mathematical 
methods are described which may be used in conjunction with the computer 
to treat linear algebraic systems. Application of these methods of treating 
linear algebraic systems to the attainment of the roots of both secular 
systems and polynomials follows. Iterative procedures for improving accu- 
racy are also included. 

A discussion of the limitations which are encountered in the REAC 
application of these methods, as well as a means of compensating for certain 
of these limitations, is included as an appendix to this paper. 

2. REAC Operational Procedures. Algebraic addition, inversion (change 
of sign), and integration with respect to time are accomplished on the REAC 
by direct current feedback amplifiers. Variables are represented by voltages; 
passage through any amplifier reverses sign. The REAC computer contains 
twenty amplifiers (7 integrating, 7 summing, and 6 inverting) and thirty 
potentiometers with the associated power supplies and ‘“‘patching’’ equip- 
ment. Each summing and integrating amplifier, without modification, can 
receive as many as seven inputs at its various input terminals. Of these 
input terminals, three have amplification factors of one, two have amplifi- 
cation factors of four, and two have amplification factors of ten. Poten- 
tiometers are used in connection with amplifier inputs-to attain other than 
the listed integral amplification factors, and, as a result, coefficients can be 
set to three place accuracy. Provision is also made for supplying initial 
condition voltages to each of the seven integrators. 

The technique frequently employed in the solution of differential equa- 
tions is as follows. First, the highest derivative is indicated as being equated 
to the other terms in the differential equation. The negative of the time 
integral of this highest derivative is anticipated at the output of a particular 
integrator. The balance of the block diagram is then set up to the point 
where the anticipation is justified in that the input to the particular inte- 
grator is the highest derivative. Figure 1 demonstrates the operation of the 
listed components; Fig. 2 presents the block diagram for the particular 
differential equation: 


dy , ,dy = 
an + 75, + 3.14y = 0, 
where at 
d 
t=0, y=10 and + = 0: 


It should further be indicated that multiplication and division of vari- 
ables are carried out on the REAC by high-speed servo units. 

The REAC was designed to serve as an electronic differential analyzer. 
We are now going to illustrate the application of such a computer to certain 
algebraic problems. 
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Fig. 1. Fig. 2. 


3. Basic Mathematical Considerations. Several methods for the ana- 
logue solution of linear algebraic systems have been investigated. Two of 
these methods are included in the discussion which follows. It will be found 
that each of these methods results in replacing the linear algebraic system 
by a linear differential system whose steady state solution exists and is the 
solution of the linear algebraic system. Obviously, in order that the two 
methods which follow be distinct, the resulting linear differential systems 
must also be somewhat different. The first of these two methods uses what 
appears to be a minimum of equipment but is usually limited in application 
to systems having a matrix of coefficients which is positive definite! [p. 110]. 
The second method to be presented is not so economical as regards com- 
puter components but has the following important advantage. The proof of 
convergence of this process of solution of linear algebraic systems is not 
dependent on the arrangement of the system of equations to be solved. 
Thus the linear system presented for solution by this method will require 
no transformations prior to being placed on the computer for solution. 

In general, given the linear algebraic system 


(1) > AyX; — B; = 0, 4=1,2,---,n 


j=1 


one can arbitrarily choose a set of X;’s such that the system can be written: 


(2) LX AgX; — Bi = &, 4=1,2,---,n, 


i=l 
where £; represents the error or residue in the ith equation. 
In the event that the matrix of coefficients is positive definite, the 
following relation has been established? [p. 22], that is, 


dX; 
a 





ana- 
ro of 
pund 
stem 
; the 

two 
tems 
what 
ation 
110]. 
com- 
of of 
s not 
Ived. 
quire 


itten: 


», the 





OTHER AIDS TO COMPUTATION 193 


or 


(3) X,=- f tae 


. The instrumentation of the following 3 X 3 system of linear equations 


x + .333y — .6672 
(4) 


8 
—.5 x+ y+.5 2 | 
25x — .25 y+ s=-7 


is set up (Fig. 3) using the relationships of equations (2) and (3). In addi- 
tion, each equation in (4) has been divided by a constant such that the 
coefficients on the main diagonal are each unity. It can be observed from 
the schematic diagram of Fig. 3 that for an mth order linear system (with 
a positive definite matrix of coefficients) the following computer components 
will, in general, be needed, 


n integrators, 
nm summing (or inverting) amplifiers, 
n* coefficient potentiometers. 
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The remaining inputs [of the total of m(m + 1) ] are each made equal to 
unity as previously indicated, thus eliminating the need for these n addi- 
tional coefficient potentiometers. 

The second method of solution of linear systems mentioned in the intro- 
ductory remarks to this section follows. It will be found that the convergence 
of this second method is not dependent on the arrangement of the system 
of equations to be solved. Hence this method is generally applicable to 
linear systems and can readily be applied to the machine methods (to be 
described later in this paper) of treating secular systems and polynomials. 

In discussing this second method for the solution of linear systems, 
attention is again focused to equations (1) and (2) of this section. These 
equations respectively define the linear system presented for solution and 
the residues resulting from choosing an arbitrary set of X,’s. The machine 
problem involved once more becomes the problem of adjusting (auto- 
matically) the arbitrarily chosen X;’s such that the residue in each equation 
will approach zero as a limit. 

It has been shown by ApcocK* and Murray? [p. 20] that if the following 
three equations are satisfied simultaneously on the computer 
5 


LAuX;-Be=t, 74 








= 1,2,---,” 
j=1 
(5) jE, = K f kde, 4=1,2,---,n¢ 
X;= - LEAs, 4=1,2,---,m” 
L ‘1 J 


each of the variables in the linear system will approach its correct value if 
a unique solution exists for the algebraic system. The existence of a unique 
solution for the algebraic system assures the existence of a steady state 
solution for the differential system (5). Due to the fact that this method is 
not restricted to the positive definite case, it is preferred by the writer and 
is referred to in the remaining sections of the paper as a general method for 
the solution of linear systems. 

It will be found that this method has been referred to in the literature 
by various descriptive titles. To this writer’s knowledge, the method has 
been referred to as the method of “minimizing the sum of the squares of 
the errors’’‘ [p. 120-121], the method of ‘‘steepest descent’’‘ [p. 45], and/or 
the ‘down the gradient adjustment process.’’* Although a numerical appli- 
cation of this process meets with some difficulty due to the indicated re- 
quirements for a means of continuous integration, this requirement is readily 
fulfilled on a differential analyzer. 

To demonstrate the routine instrumentation of this second method of 
solution, the following third order linear system has been chosen: 


0.750X, — 0.500X2 + X; — 32.55 =0 
0.400X; + X2 — 0.600X; + 18.0 =O. 
Xi + 0.500X2 + 0.333X; — 11.66 = 0. 


For convenience in setting the system on the computer, each equation has 
been divided by a constant such that the magnitude of the largest coefficient 


(6) 
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in each equation is < 1. Following this procedure, it becomes necessary to 
satisfy simultaneously the relations (5). The schematic diagram for achiev- 
ing the solution is in Fig. 4. 

The first approximations to the solutions attained by the method of 
Fig. 4 were found to be accurate to two significant figures. It should be 
pointed out, however, that accuracy considerations are a function of the 
system of equations being solved as well as of the precision of the compo- 
nents being used. In general, unless the equations are so ill-conditioned* 
[p. 45] that the approximations made in the setting of coefficients affect 
the solution considerably, a solution to any required accuracy can in prin- 
ciple be obtained by use of the iterative procedures to be described in 
Section 4. 

















Fig. 4. 


An inspection of Fig. 4 makes evident the fact, as does the system (5), 
that this method requires a duplication of the coefficient network. Murray 
has indicated that if an adjusting type of machine is to operate successfully 
whenever the determinant A is not zero, the feedback system will of neces- 
sity involve the use of the A,; twice.5 

Further, a comparison of the schematic diagrams, Fig. 3 and Fig. 4, 
demonstrates that the general method for soiution of linear systems, Fig. 4, 
requires approximately twice the machine components as does the method 
for the positive definite case, Fig. 3. This increase in component require- 
ments of the general method leads to certain limitations of the method. 
These limitations, as well as a suggested (and demonstrated) method of 
alleviating certain of the inherent difficulties, are to be found in the appendix. 
4. Iterative Procedures for Extending Accuracy. In the event that the 
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computer-attained solutions for a linear system are not as precise as desired, 
the following iterative process* may be used to extend this accuracy. These 
initial REAC-attained values of the roots may be substituted into the origi- 
nal set of equations (on standard desk calculator) and the residues, R;, 
computed for each equation. These R,’s (¢ = 1,2, ---,) are then each 
multiplied by a constant, the magnitude of which is chosen to utilize as 
much scale factor as possible. These multiples of the residues are used as 
new constants in the system of equations (for example, 100 R; replaces b; 
(4 = 1,2, ---, )) and the system is again solved on the REAC. This process 
requires only that the coefficients on the computer representing the ),’s 
(4 = 1,2, ---, 2) be reset in order that these increments to the roots can be 
attained. These new values, (AX;,’s), are then combined with the initial set 
of REAC-attained roots to yield a more precise set of roots. This iteration 
process can be repeated if still higher accuracy is desired. 

5. Linear Algebraic Systems with Complex Coefficients. Linear alge- 
braic systems with complex coefficients such as those arising in certain 
aeroelastic investigations may also be solved by this general method. The 
process which has been used doubles the order of the system by separating 
the equations into their real and imaginary parts. The pairing of the real 
and imaginary parts of the solution is provided for as a result of the relation 
that 2; = x; + 14; G = 1,2, ---, 2). 

6. Eigenvalue Problems. Having available a machine method (the gen- 
eral method, section 3) for solving a set of linear algebraic equations which 
has a unique solution, attention is given to the application of this general 
method to the solution of eigenvalue problems. It will be shown that an 
adaptation of this method can be used for the set of equations, 


(7) D> (Ag — djA)X; = 0, $= Te & oe) Ts 

7=1 
In order that the system (7) have solutions other than the trivial 
X, = X, = --- = X, = 0, it is necessary and sufficient that the determi- 
nant of the coefficients vanishes, that is, 


(Au — byd)(Aiz — Dd) «++ (Ain — Bind) 
(8) (An — Dud) (Ars - baad) - “An i Dank) 


(doa _ b ar) (Ans le b od) - ven — Bank) 


Now, besides the pre-assigned quantities A;; and 5;;, equation (8) involves 
only the parameter A, so that (7) has non-trivial solutions if and only if \ 
satisfies equation (8). Equation (8)—an equation of the mth degree in the 
parameter \—is called the secular equation of the matrix,’ and we have the 
result: the eigenvalues of a matrix are the roots of the secular equation of 
the matrix. Although in the numerical process the eigenvalues of a given 
matrix can be determined directly through the secular equation (8), and 
without reference to the defining equation (7), this method is not followed 
in the machine set-up presented here. 

For the purpose of simplifying the presentation of the recommended 
machine method, a 3 X 3 system having all real roots will be used. It is a 
well known fact that the eigenvalues of a symmetrical characteristic equa- 
tion are all real. Hence, the quotation from GurLLemrn® that [p. 112] ‘‘the 
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matrices encountered in practical problems are predominantly symmetrical” 
adds significance to a machine method which can be used for the determina- 
tion of real roots. This is true even though the determination of real roots 
is a special case of the more general problem of determination of both real 
and complex roots. This general problem will be discussed following the 
presentation of the machine method for real roots. 

The eigenvalue problem (having all real roots) which will be used as an 
example follows: 


{- (.611 — A)X, — .333X2 + .111X; 


(9) .333X, + (.556 — A)X2 — .222X; 


0 
0 
-111X, — .222X, + (.333 — r»)X; = 0 


This system (9) has non-trivial solutions if, and only if, \ assumes a proper 
value. If \ assumes a proper value, and if the roots of the system are all 
distinct, the rank of the determinant (10) 


(.611 — A) — .333 ll 
(10) — .333 (.556 — A) — .222 
111 —.222 (.333 — d) 








will be m — 1 (in this case m — 1 = 2). Hence one of the variables X; will 
be arbitrary, and it is this fact which becomes the basis for the machine 
method which follows. 

Thus the first step in setting up the system for REAC solution results 
in assigning to one of the variables X; of (9) an arbitrary value. In this 
example X; is assigned the value 10, and the system becomes 


(.611 — A)X, — .333X; + 1.11 = 0 
(11) — .333X, + (.556 — A)X_ — 2.22 = 0 


-111X, — .222X2 + (3.33 — 10A) = 0 


Following this procedure, the last two equations of (11) may be set up for 
automatic solution on the computer with provision made for varying A 
(Fig. 5). 

For any given value of \ the computer automatically solves these last 
two equations for X; and X_. The terms in the first equation are continu- 
ously summed in a separate amplifier (amplifier number 14) whose output 
indicates the error or residue of the first equation. Hence, as A is (slowly) 
varied, this error term of the first equation may be observed by the operator. 
When for a given setting of \ this error becomes and remains zero, a “‘proper”’ 
value or eigenvalue has been determined. 

Attention is called to the fact that in the schematic (Fig. 5), the simul- 
taneous variation of \ throughout the system was accomplished by utilizing 
REAC multiplier servos. There would be advantages to incorporating man- 
ually operated ganged potentiometers in place of using these servos. How- 
ever, ganged potentiometers were not available when these machine methods 
were tested. 

7. Iterative Procedures for Improving Accuracy. Iterative procedures 
are available for improving the accuracy of the eigenvalues and eigen- 
vectors.® If \! and X! are approximations to a solution of 


(12) (A — \B)X = 0, 





(17) (on — Mbu)AX, it+(n - Mba) AXa+ -- "+204 = ~ Aiba, n1)AX, n-1+0 
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Fig. 5. 


an equation of the same general form as (7), then 

(13) [A — (A! + AN)BILX! + AX] = 0. 

This equation (13) upon expanding becomes 

(14) [A — \'B]X!+ [A — \BJAX — AABX! — ANBAX = 0. 

Now if (A — \'B)X' is accurately computed (on a desk calculator) it then 
becomes the matrix of error 

(15) [A — NB]X! = —e 

and furthermore, if the second order term AX BAX is neglected, equation (14) 
can be written 

(16) [A — NB]JAX — A\BX! = «. 

As was previously pointed out, in the discussion of the computer method 
for handling this type of problem, X, is considered to be arbitrary and is 


set equal to any convenient constant. As a result, AX, is taken to be zero 
and equation (16) can be expressed as 


(Qu — MO) AXi+ (Giz — MD) AXo+ + +» + (Gina — MM, na AXn st “| 41 
—Ad) Ke 


= sl 


as —d! bu) AX, + (ans - Nbna)AXa+ «+ "e den a-1— ~ Ntbp, a ae +0 


Ka en 
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where 


K; = baX,! +- bX! +--+ + binX aI. 
Equation (17) when written in the detailed form of a simultaneous linear 
algebraic system becomes 
(dix — XB) AX4+ (Giz — AD) AX2+ + + + + (G1, na — Mi, na) AXn a — KiAA= 
(18) (an ms Mbn)OX1+ (a - Nba)AXs+ - - “+(a ann = dibs, na)AXs a —Ksi= es 


es —Nbas)AXs + (Gna —NMbgs)AX+- « ‘Silas 2-1 Abs, a ~K.Ad =€,n 


Hence, this ordinary linear system (18) may be solved on the analogue 
computer for the unknown increments AX,, AX», ---, AX,-1, and AX. It is 
noted that after the e; and K; have been computed (on a desk calculator), 
only the following two steps are necessary prior to actuating the computer 
to attain the desired solution yielding the listed incremental values. 


(1) The last column of coefficients (on the computer) is to be readjusted 
to the K; values and 

(2) The column of constants (originally set to zero) is to be adjusted to 
the values ce; where c regulates the scale factor. 


The approximation \! remains on the computer as a constant and the 
computer automatically determines the increment Ad as well as the incre- 
ments to each of the variables X; (¢ = 1, 2, ---,m — 1). 

8. Considerations for Complex Eigenvalues. The procedures, previously 
outlined for determining real eigenvalues and eigenvectors, can be extended 
to the determination of the complex eigenvalues and eigenvectors of sys- 
tems having either real or complex coefficients. However the method of 
treatment effectively halves the capacity of the machine. This results from 
the fact that the manner of handling the system, in substance, replaces an 
n variable complex system by a 2m variable real system. The procedure to 
be followed after having converted to this real system follows the basic 
method used in setting up Fig. 5. However, the following two details must 
be kept in mind. The arbitrary variable will be one of the original complex 
variables, €.g., Zn = Xn + tyn, and may be assigned a convenient arbitrary 
value (such as z, = 1). Secondly, the parameter A is also complex and of 
the form \ = u + iv. Hence, provision must be made for independently 
varying both the real and imaginary parts of the parameter A of the original 
complex system. 

9. Application of the Method of Secular Systems to Polynomials. In 
view of the previous discussion (sections VI and VIII), it is evident that 
the roots of polynomials could be found by an application of these methods 
if the polynomial can be written in the form of (7). This process can be 
accomplished and the details of the method are exhibited for the following 
fourth degree polynomial. 

Given the polynomial 


At + DA? + BA? + BA +, = O 
it, in turn, can be written in the following form 


inet oF be 
a -% Ss bs | _ 
(19) ee b, | = 2 


0 0 -1 (A+)d) 
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By comparison, one can see that this is a special case of equation (8). 
To make possible the machine set-up of this equation (19), we utilize a set 
of auxiliary variables ¢,, t2, ts and #& in a manner such that the resulting 
system becomes 


ty a. bg ts = 0 

— ty +Abe + bts = 0 

(20) —™ to + Atz + be tg = 0 
— t+(A+d)te = 0 


This system of equations (20) can be seen to be a special case of the more 
general system (7). Hence the procedures described for eigenvalue problems 
can be utilized in setting this system up for machine solution. This being 
the case, the method can be used to attain the roots of polynomials depend- 
ing, of course, on the degree of the polynomial. 

10. Appendix. It was pointed out in this paper that a REAC utilization 
of the general method (section 3) meets with limitations. However, we can 
compensate for certain of these limitations. The following paragraphs indi- 
cate one method whereby it seems feasible to use a pair of REAC computers 
in attaining the solution of a 10 X 10 linear algebraic system with real 
coefficients. It also seems likely that one could surmount the difficulties 
involved in extending these measures to the utilization of three REAC’s in 
solving a 15 X 15 linear system. In such,an event, the ‘“‘condition” of the 
system of equations becomes more important as regards accuracy consid- 
erations. This would be important in view of the precision of the components 
to be used. 

There are two basic problems involved in the instrumentation of this 
method in conjunction with the REAC. The first problem pertains to a 
satisfactory means of providing a sufficient number of coefficient poten- 
tiometers. Allied with this is the necessity of preventing amplifiers from 
‘limiting’? as a result of the number of loads an amplifier may be called 
upon to supply. The second problem is to determine a suitable method of 
patching the problem to the computer. For a method to be suitable, it must 
reduce not only the time and effort involved in patching, but it must also 
reduce the likelihood of human error in the performance of these tasks. 

The method of handling the first problem mentioned above includes the 
following considerations. The output impedance of a REAC amplifier is 
about 5KQ. In general, each amplifier has to feed a number of coefficient 
potentiometers; for a given amplifier, this number will be equal to the 
number of coefficients of the same sign in a given row (or column) of the 
matrix of coefficients. Since the maximum of this number is ten, it is neces- 
sary to use potentiometers of rather high resistance, so that ten such 
potentiometers in parallel will not overload an amplifier. Suitable wire 
wound potentiometers are not available at a reasonable price and, as a 
result, a combination of a resistor decade and a potentiometer" is used to 
provide the equivalent of a coefficient potentiometer having approximately 
a 180KQ total resistance (Figs. 6 and 7). Each of these potentiometer 
combinations is set to its proper ratio by applying a constant voltage to 
the input and then measuring the voltage at the output of an amplifier 
which is fed from the potentiometer arm. Thus the need for calibrated 
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components (precision resistors and calibrated linear potentiometers) in 
conjunction with a load correction table is eliminated. 

As a result, it is believed that an adaptation of the above-mentioned 
resistor decade and potentiometer combination provides a means of satis- 
fying the impedance requirements of an amplifier when it may be necessary 
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Fig. 6. 


for that amplifier to supply as many as ten output loads. Further, this 
method of furnishing voltage dividers to be used for the setting of coeffi- 
cients is found to be, comparatively, quite economical. 

If such a res'stor decade and potentiometer combination is accepted as 
a means of alleviating the first problem of instrumentation cited above, one 
may turn to the second problem stated. This second problem regards itself 
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with the “patching” of the problem to the computer. The suggestions 
followed in this case are based on the knowledge that the general method 
is stable and convergent for all linear systems that are not too “‘ill-condi- 
tioned.’’ Hence, with the exception of setting the proper numerical values 
of the coefficients, the wiring diagram (the problem of “‘patching’’) is 
identical for all linear systems of a given order. This fact, along with the 
intention of providing for systems up to and including a 10 X 10 linear 
algebraic system, led to the following type instrumentation. The devices 
previously described for the setting of coefficients are arranged in the form 
of the matrix coefficients of a linear system. This arrangement is mounted 
(Fig. 9) such that it can easily be moved near to or away from the REAC 
installation. All wiring (see schematic Fig. 8) of the matrix network is 
incorporated by permanent connections (soldered leads) with the exception 
of the leads from the buss lines to the respective grid terminals and output 
terminals of the REAC amplifiers. With the present availability of remov- 
able ‘‘patch bays’’ for the REAC, the connections from the buss lines to a 
set of patch bays can also be made permanent if desired. 

In such a case, after the coefficients have all been properly adjusted, 
the attachment of these patch bays to the computer is the only requirement 
prior to actuating the computer to attain a solution. Hence, the REAC 
installation is required only for the time necessary for the method to con- 
verge (with ¢ = time in seconds) along with the time required to read the 
values of the variables. Iterations for improving accuracy, if necessary, will 
require additional REAC time. 





Fig.7. 
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Several 10 X 10 linear algebraic systems were solved on the REAC 
installation using the above-described coefficient matrix panels in conjunc- 
tion with the general method of section 3. These examples indicated that 
the full accuracy of the REAC may be realized in satisfying simultaneously 
the individual equations of the linear system. The accuracy of the solutions 
is, of course, dependent upon the ‘‘condition’’ of the system of equations 
being solved. 


The author wishes to express his thanks to Dr. M. G. SCHERBERG, Chief, 
Applied Mathematics Research Team, for initiating and making possible 
this investigation and to Dr. E. P. Littte, Mr. L. M. WARSHAWSKY, and 
Mr. W. Braun who also assisted with their criticism and advice. 


LANDIS GEPHART 
Flight Research Laboratory 
Wright-Patterson Airforce Base 
Ohio 


1D. V. Wipper, Advanced Calculus. New York, 1947. 
=e F. J. Murray, The Theory of Mathematical Machines, revised edition. New York, 

7W. A. ADCOCK, “Automatic simultaneous equation computer and its use in solving 
secular equations,” Rev. Sci. Instruments, v. 19, 1948, p. 181-187. 

*D. R. HARTREE, Calculating Instruments and Machines. Urbana, 1949. 

5 F. J. Murray, “Linear equation solvers,” Quart. Appl. Math., v. 7, 1949, p. 263-274. 

SE. A. GOLDBERG, RCA Rev., v. 9, 1948, p. 394-405. 

7V. ROJANSKY, Introductory Quantum testa. New York, 1947. 

8 E. A. GuILLEMIN, The Mathematics of Circuit Analysis. New York, 1949. 

*P. A. Dennis & D. G. DILL, ‘Application of simultaneous-equation machines to 
aircraft structure and flutter problems,” Jn. Aero. Sci., v. 17, 1950, p. 107-11 

M. WarsHAwsky & W. Braun of the Flight Research Ft a WPAFB, 

ret the design of this resistor decade and potentiometer combination as well as the 
design of the coefficient matrix panels pictured in the frontispiece. 





204 NOTES 


NOTES 


136. SOLUTION OF y = e’.—In Note 20 [MTAC, v. 1, p. 202-203] a 
quotation from the Mathematical Gazette regarding a 60D approximation to 
the smaller solution of the equation 10 log y = y, or y = 10¥/! = 1.259y, 
was the basis of varied comments. Note 25 [MTAC, v. 1, p. 334-335] 
showed that EuLER had given a 7D approximation to this same y. In 
Note 38 [MTAC, v. 1, p. 431 ] a quotation from an early volume of Institute 
of Actuaries, Journal, was given without comment as a curiosity, since it 
suggested that a solution of y = e” might be found by computation from 
the series 
4-4-4 5-5-5 
2-3-4 + 2-3-4°5 





3°3 
yelt+; “+ 2.3 + 
That this series did not converge, the author of the quotation did not 
mention. 

For six years since Note 38 was published, no reader wrote to me re- 
garding the obvious error in this quotation. But finally, Mr. Cec. HAsTINGs, 
Jr., of Santa Monica, Calif., has now incidentally commented upon it in a 
letter. 

In RMT 884 [MTAC, v. 5, p. 140] tables and graphs of the equation 
y = x” were considered, and it was clear that there were real values of y 
only if x < e!/* = 1.44, and therefore in particular for x ~ 1.259. Hence, 
for x = e, y is necessarily transcendentally imaginary. (Analogous to the 
equation in Note 20 we have the equation eln y = y, for which the two 
solutions coincide in y = e.) 

R. C. ARCHIBALD 
Brown University 
Providence, R. I. 


137. PUNCHED CarD TABLEs.—In order to keep our Guide to tables on 
Punched Cards [MTAC, v. 5, p. 185-212] up to date, it is planned to 
collect new material and errata for insertion in a future number of MTAC. 
Readers are requested to send to the undersigned any information relative 
to punched card tables. 

Several laboratories on the West Coast are collaborating in the punching 
of basic tables during spare time. In many laboratories which do accounting 
work in addition to mathematical computing, there are peak periods in the 
month when all key-punch operators are busy, and these periods are followed 
by one or two weeks of relatively light work. During these latter periods 
key punch time could be used to prepare important card tables from the 
British Association for the Advancement of Science tables, and the older 
ones of the Mathematical Tables Project which are not available on punched 
cards. It might be mentioned that the Table of Arc sin x was punched on 
this basis by BENJAMIN FERBER of Consolidated Vultee Aircraft Corpora- 
tion, and that he also found time to difference the key-punched values. 
Similarly, CHARLES Davis of the North American Aviation Corporation 
contributed the Tables of Circular and Hyperbolic Tangents and Cotan- 
gents, which he checked by multiplication of the reciprocal functions. If 
other laboratories can undertake a similar service, the undersigned are 
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prepared to act as a “clearing house,” supplying information about the 
most urgent tasks that call for key-punching. 

G. BLANCH 
E. C. YOWELL 
NBSINA 


138. A NEW MERSENNE PRIME.—The program described in Note 131 (c) 
has produced the 15th Mersenne prime: 2'* — 1 on June 25. The SWAC 
tests this number in 13 minutes 25 seconds. 


D. H. L. 
NBSINA 


QUERY 


42. Where may one obtain a set of the BEEVERS and Lipson strips for 
the summation of Fourier series? These are described in the following 
references: C. A. BEEVERS & H. Lipson, “A rapid method for the summa- 
tion of a two-dimensional Fourier series,’ Phil. Mag., s. 7, v. 17, 1934, p. 
855-859, G. L. CLARK, Applied X-Rays. New York, 1940, p. 322. 

J. M. WaITE 


Magnolia Petroleum Co. 
Box 900 
Dallas, Texas 


QUERIES—REPLIES 


49. TABLE OF MULTIPLICATION (Q. 40, v. 6, p. 61).—I possess the first 
edition of the tables of J. B. Ovon, bound in one volume, giving the products 
of all the integers up to 509 X 500 without mention of the author’s name. 
The title is: 

Tables de multiplication, 4 l'usage de MM. les ingénieurs employés au 
cadastre de la France, et de MM. les directeurs des contributions, chargés de la 
vértfication des opérations arithmétiques des ingénteurs. 

Approuvées par S.E. le Ministre des Finances 

(Prix: 15 francs) 


A PARIS 
a l’imprimerie de VALADE, rue Coquilliére 


AN XIII (1805) 


The first page bears the following annotation written in pencil: Renard, 
rue Ste Anne 71, the second tome 15 f by M. Oyon. 

I have no information about Oyon (he may have been an engineer in 
the Land Registry Service). 


R. LréNarD 
Tulle (Corréze) France 
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